JP3882584B2 - 乱流の解析システム並びに乱流解析プログラム - Google Patents
乱流の解析システム並びに乱流解析プログラム Download PDFInfo
- Publication number
- JP3882584B2 JP3882584B2 JP2001340119A JP2001340119A JP3882584B2 JP 3882584 B2 JP3882584 B2 JP 3882584B2 JP 2001340119 A JP2001340119 A JP 2001340119A JP 2001340119 A JP2001340119 A JP 2001340119A JP 3882584 B2 JP3882584 B2 JP 3882584B2
- Authority
- JP
- Japan
- Prior art keywords
- component
- scale
- model
- sgs
- stress
- 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.)
- Expired - Fee Related
Links
Images
Description
【０００１】
【発明の属する技術分野】
本発明は流体機器や流体関連機器における流れの数値解析方法に関わり、特に乱流場を非定常解析する計算方法であるラージ・エディ・シミュレーションにおいて、ＳＧＳ応力のモデル係数を自動的に決定する方法に関する。
【０００２】
【従来の技術】
流体機器の設計や開発では、流れの数値解析が利用されることが多い。そして対象とする流れ場の多くは乱流である。乱流の解析方法の一つにラージ・エディ・シミュレーション（ＬＥＳ）がある。
【０００３】
流れの数値解析手法の多くは、解析の対象となる領域を、空間的に有限個の領域、あるいは有限個の物理量定義点に分割して扱う。これを空間の離散化といい、一般的にはメッシュ分割，要素分割，セル分割などとよぶ。空間を離散化すると、たとえばメッシュ幅よりもサイズが小さい渦などの流れ構造は、もはや解像できなくなってしまう。しかし、こうした微小な流れ構造も、全体の流れ場に影響を及ぼすので、その影響を考慮する必要がある。こうした微小な流れ構造（波数空間上での高波数域の速度の変動成分）の影響をサブグリッドスケールモデル（ＳＧＳモデル）として考慮し、流れを非定常計算するのがＬＥＳである。
【０００４】
空間の離散化は、物理量の空間的な変動の、高波数成分を除去する一種のローパスフィルタとみなすことができ、このフィルタをグリッドフィルタとよぶ。波数空間において、グリッドフィルタにより除去される高波数の領域のことをサブグリッドスケールとよび、除去されない低波数域のことをグリッドスケールとよぶ。また波数空間上でグリッドスケールに含まれる物理量の空間変動成分のことをグリッドスケール成分（ＧＳ成分）とよび、サブグリッドスケールに含まれる空間変動成分をサブグリッドスケール成分（ＳＧＳ成分）とよぶ。
【０００５】
任意の物理量をφであらわすとき、空間的に分布するφにグリッドフィルタを作用させて得られる値を
【０００６】
【数４２】
【０００７】
であらわす。
【０００８】
空間の離散化で解像できなくなった微小な流れ構造、すなわち空間的な速度変動のＳＧＳ成分が、流れの運動量変化に及ぼす影響は、応力の形であらわれる。この応力をサブグリッドスケール応力（ＳＧＳ応力）という。
【０００９】
速度のｉ方向成分をｕi と表すとき、ＳＧＳ応力は数１であらわされる。
【００１０】
【数１】
【００１１】
ＳＧＳ応力中の
【００１２】
【数４３】
【００１３】
については直接的には計算することができないため、ＳＧＳ応力に対してモデル化が必要となる。
【００１４】
ＳＧＳ応力に対する代表的なＳＧＳモデルに、スマゴリンスキーモデルがある。スマゴリンスキーモデルは数２で表される渦粘性型のＳＧＳモデルである。
【００１５】
【数２】
【００１６】
ここでグリッドスケールでの歪率テンソルは数３で与えられる。
【００１７】
【数３】
【００１８】
スマゴリンスキーモデルでは渦動粘性係数をモデル定数Ｃを含んだ数４で与える。
【００１９】
【数４】
【００２０】
なお、本明細書においては数式の記述において、アインシュタインの総和規約を用いるものとする。また任意のテンソルΦijに対し、Φij * はその非等方成分Φij *≡Φij−δijΦkk／３をあらわすものとする。
【００２１】
スマゴリンスキーモデルのモデル定数は流れの性状によりその最適値が異なる。例えば一様等方性乱流ではその最適値は０.０３ から０.０４ 程度の値であり、チャンネル乱流では０.０１ 程度の値が最適値である。
【００２２】
そのためスマゴリンスキーモデルを用いた解析で精度の高い結果を得るためには、対象とする流れ場に適したモデル定数を用いる必要がある。しかし複雑な流れ場を対象とするとき、そのようなモデル定数のチューニングは極めて難しい。
【００２３】
そこで流動状態に応じてモデル定数の局所的な値を自動的に算出する方法が提案されている。その方法はフィジックス オブ フルイズ エー ４（１９９２年）第６３３頁から第６３５頁（Physics of Fluids A 4 (1992) pp.633-635 ）に記載されたものである。この方法をもちいたＳＧＳモデルはダイナミック・スマゴリンスキーモデル（ＤＳＭ）と呼ばれる。またこのときモデル定数を算出するために用いる計算手続きはリリーの最小自乗法とよばれる。
【００２４】
ＳＧＳモデルにスマゴリンスキーモデル以外のモデルを適用する場合についても、モデル定数を自動的に算出する方法が提案されている。例えば数５はスケール相似則モデルの一つであるバーディナモデルであるが、これと数２，数４のスマゴリンスキーモデルとを線形結合したモデルとして、数６に示すミックスドモデルがある。
【００２５】
【数５】
【００２６】
【数６】
【００２７】
数６のミックスドモデルを対象として、式中に含まれるスマゴリンスキーモデル部のモデル定数Ｃを算出する方法が、フィジックス オブ フルイズ Ａ ５（１９９３年）第３１８６頁から第３１９６頁（Physics of Fluids A 5(1993)pp.3186-3196）に記載されている。この方法はダイナミック・ミックスドモデル（ＤＭＭ）とよばれる。ダイナミック・ミックスドモデルにおいても、モデル定数の算出手続きにはリリーの最小自乗法が用いられている。
【００２８】
また数７に示すツーパラメータ・ミックスドモデルはミックスドモデルの一種であり、モデル定数を２つ含むモデルである。
【００２９】
【数７】
【００３０】
この２つのモデル定数を自動的に算出する方法はフィジックス オブ フルイズ ７（１９９５年）第２８３１頁から第２８４７頁（Physics of Fluids 7 (1997) pp.2831-2847)に記載されている。この方法はダイナミック・ツーパラメータ・ミックスドモデル（ＤＴＭＭ）と呼ばれる。またこの方法でモデル定数を算出するために用いられる方法はやはり、リリーの最小自乗法である。
【００３１】
【発明が解決しようとする課題】
上記の従来技術における問題点は、数値的に不安定な点である。すなわち、ＳＧＳモデルのモデル定数の自動算出に従来の技術を用いた解析では、数値振動が生じやすく、場合によっては非物理的な振動解が計算結果として得られる。さらには計算途中においてゼロ割りやフローティングポイントオーバーフローなどが生じ、計算の続行が不能となることもある。
【００３２】
こうした数値不安定を防ぐ目的で、統計的に一様な方向への平均操作を用いる方法が提案されており、その一例はフィジックス オブ フルイズ エー ３（１９９１年）第１７６０頁から第１７６５頁（Physics of Fluids A 3 (1991)pp.1760-1765）に記載されている。但しこうした平均操作は単純な体系にしか適用できず、またその物理的な根拠も不明であるため、一般的な解決策ではない。乱流の解析方法の一つであるＬＥＳにおいて、ＳＧＳ応力のモデル式中に含まれるモデル係数もしくはモデル定数の局所的な値を、瞬時の流動状況に応じて自動的に算出することが望ましい。さらに加えて、数値的な不安定現象が生じにくい、数値的に安定な計算方法であることが望ましい。
【００３３】
本発明の目的は、乱流の解析方法の一つであるＬＥＳにおいて、ＳＧＳ応力のモデル式中に含まれるモデル係数もしくはモデル定数を、自動的に算出し、数値的に安定な乱流の解析方法，乱流の解析システム並びに乱流解析プログラムを提供することにある。
【００３４】
【課題を解決するための手段】
本発明は、該離散化に伴い空間的には解像することができなくなった高波数域の速度の空間変動成分が流体の運動に与える影響であるサブグリッドスケール応力（ＳＧＳ応力）をモデル化したモデル式を計算する計算手順が前記計算機メモリもしくは記録媒体に保存され、物理量の空間変動の波数成分のうち該離散化に伴い空間的には解像することができなくなった高波数域 ( サブグリッドスケール ) の成分をサブグリッドスケール成分（ＳＧＳ成分），物理量の空間変動の波数成分のうちＳＧＳ成分ではない低波数域（グリッドスケール）の成分をグリッドスケール成分（ＧＳ成分）とし、物理量の空間変動の高波数成分を除去する作用をもつ空間フィルタであるテストフィルタを、物理量の空間変動の前記
ＧＳ成分に対して作用させ、前記ＧＳ成分のうちの高波数域の成分を除去した残りの低波数域（テストスケール）の成分をテストスケール成分（ＴＳ成分）とし、物理量の空間変動の波数成分のうち前記ＴＳ成分には含まれない高波数域（サブテストスケール）の成分をサブテストスケール成分（ＳＴＳ成分）とする場合、速度ベクトルの前記ＧＳ成分を算出し、前記計算メモリもしくは前記記録媒体に記憶する手段と、前記算出された前記速度ベクトルのＧＳ成分のデータに基づき、運動エネルギーのＴＳ成分の粘性散逸量の近似値と、運動エネルギーのＧＳ成分の粘性散逸量の近似値と、テストスケールからサブテストスケールへのエネルギー伝達量の近似値とを算出する手段と、運動エネルギーの該ＴＳ成分の粘性散逸量とテストスケールからサブテストスケールへのエネルギー伝達量との和が、運動エネルギーの該ＧＳ成分の粘性散逸量とグリッドスケールからサブグリッドスケールへのエネルギー伝達量との和に等しくなるという、前記計算メモリもしくは前記記録媒体に保存された局所スケール間平衡仮定関係式の計算手順を呼び出す手段と、グリッドスケールからサブグリッドスケールへのエネルギー伝達量を、ＳＧＳ応力テンソルと、速度の空間変動のＧＳ成分に対する歪率テンソルとの内積に（−１）を乗じたのち、テストフィルタを作用させて得られる値であるとみなし、このとき該ＳＧＳ応力テンソルを、その値をまだ定めていないモデル係数を、１つだけ含む前記モデル式で与えるとき、該局所スケール間平衡仮定関係式を、該値をまだ定めていないモデル係数について整理することにより、該ＳＧＳ応力のモデル式に含まれる該モデル係数を算出するためのモデル係数算出式を導き、該モデル係数算出式を用いて該モデル係数を計算してグリッドスケールからサブグリッドスケールへのエネルギー伝達量を計算する手段とを備えたことを特徴とする。
【００３５】
【発明の実施の形態】
ＳＧＳ応力は流れの運動量ベクトルを時間・空間的に変化させる作用をもつ。そのため必然的に、流れの運動エネルギー変化に対しても作用を及ぼすことになる。
【００３６】
数８にサブグリッドスケール運動エネルギー（ＳＧＳ運動エネルギー）の時間変化に関する簡略式を示す。
【００３７】
【数８】
【００３８】
ＳＧＳ運動エネルギーは、速度の空間変動のＳＧＳ成分がもつ運動エネルギーに相当する。数８の中の保存項は、移流項や拡散項など、ＳＧＳ運動エネルギーの輸送に関連した項であり、ＳＧＳ運動エネルギーの生成や消滅には直接的には関係しない項である。ＳＧＳ運動エネルギーの生成・消滅に関与するのは、運動エネルギーのＳＧＳ成分の粘性散逸量と、グリッドスケールからサブグリッドスケールへのエネルギー伝達量である。数８において、粘性散逸は消滅項として働き、グリッドスケールからサブグリッドスケールへのエネルギー伝達は生成項として働く。運動エネルギーのＳＧＳ成分の粘性散逸量は数９で与えられ、運動エネルギーのＧＳ成分の粘性散逸量は数１０で与えられる。
【００３９】
【数９】
【００４０】
【数１０】
【００４１】
また、グリッドスケールからサブグリッドスケールへのエネルギー伝達量は数１１で与えれれる。
【００４２】
【数１１】
【００４３】
数１１に示されるように、ＳＧＳ応力が関与するのは、グリッドスケールからサブグリッドスケールへのエネルギー伝達量に関してである。
【００４４】
数１２にグリッドスケール運動エネルギー（ＧＳ運動エネルギー）の時間変化に関する簡略式を示す。
【００４５】
【数１２】
【００４６】
ＧＳ運動エネルギーは、速度の空間変動のＧＳ成分がもつ運動エネルギーに相当する。数１２において、運動エネルギーのＧＳ成分の粘性散逸は消滅項として働き、グリッドスケールからサブグリッドスケールへのエネルギー伝達量もまた消滅項として働く。すなわち、ＧＳ運動エネルギーの一部が消滅し、その消滅分がＳＧＳ運動エネルギーの生成分となる。
【００４７】
図１に、上に述べた波数空間での運動エネルギー伝達に関する概念図を示す。図１は波数空間上でのエネルギースペクトルを表しており、運動エネルギーの
ＧＳ成分の粘性散逸１、およびＳＧＳ成分の粘性散逸２により、運動エネルギーの一部は熱として散逸する。また、運動エネルギーのＧＳ成分の一部は、グリッドスケールからサブグリッドスケールへのエネルギー伝達３により、運動エネルギーのＳＧＳ成分に移動する。
【００４８】
上述した波数空間上での運動エネルギーの伝達特性は、計算の数値安定性および計算精度の観点から極めて重要である。グリッドスケールからサブグリッドスケールへのエネルギー伝達量を過小評価するということは、ＧＳ運動エネルギーの減衰が小さいということで、グリッドスケールの高波数域に運動エネルギーが蓄積されることになる。すなわち、非物理的な高波数の数値振動が生じやすくなる。その反対にエネルギー伝達量を過大評価すると、ＧＳ運動エネルギーの減衰が大きくなり、特に高波数域での減衰が大きくなる。そのため微小な流れ構造を解像することができなくなってしまい、解析結果の精度が悪くなる。したがってＳＧＳ応力のモデルには、グリッドスケールからサブグリッドスケールへのエネルギー伝達量を適切に模擬することが求められる。
【００４９】
従来の技術に記載したＳＧＳモデルのモデル定数の決定法はすべて、リリーの最小自乗法を用いたものであった。ここでリリーの最小自乗法について説明する。
【００５０】
まず、物理量の空間変動の高波数成分を除去する空間フィルタであるテストフィルタを導入する。このとき任意の物理量φにテストフィルタを作用させて得られる値を
【００５１】
【数４４】
【００５２】
と表記することにする。また、テストフィルタを物理量のＧＳ成分に作用させ、ＧＳ成分のうちの高波数域の成分を除去した、残りの低波数領域をテストスケールとよび、その成分をテストスケール成分（ＴＳ成分）とよぶ。また全波数域からテストスケールを除去した残りの高波数領域をサブテストスケールとよび、その成分をサブテストスケール成分（ＳＴＳ成分）とよぶことにする。
【００５３】
このとき、空間的な速度変動のＳＴＳ成分が運動量に及ぼす影響であるサブテストスケール応力（ＳＴＳ応力）は、数１３で与えられる。
【００５４】
【数１３】
【００５５】
ＳＧＳ応力とＳＴＳ応力の間には、数１４に示す関係式が成り立つ。
【００５６】
【数１４】
【００５７】
この関係式は恒等式であり、ジャーマノ恒等式と呼ばれる。
【００５８】
リリーの最小自乗法では、ジャーマノ恒等式中のＳＧＳ応力およびＳＴＳ応力に対し、同一のＳＧＳモデルを適用する。また、このときの残差テンソルを数１５で定義する。
【００５９】
【数１５】
【００６０】
ＳＧＳ応力に対するモデル式にＣ1からＣnまでのｎ個のモデル定数が含まれる場合について考える。リリーの最小自乗法では、数１６に示すように、残差テンソル成分の自乗和の、モデル定数に対する偏微分値が、すべてのモデル定数に対して０となるよう、モデル定数Ｃ1〜Ｃnを求める。
【００６１】
【数１６】
【００６２】
この手続きは一般の最小自乗法と同じである。すなわち、リリーの最小自乗法は、ジャーマノ恒等式をなるべく満足するように、モデル式で与えられたＳＧＳ応力テンソルの成分値を最適化する方法である。
【００６３】
このように、リリーの最小自乗法自体は、グリッドスケールからサブグリッドスケールへのエネルギー伝達特性に関しては、直接的には何一つ言及していない。しかし、ジャーマノ恒等式に関する残差テンソルが十分に小さくなるならば、ＳＧＳ応力テンソルの成分に関する近似精度が十分に高くなることが期待でき、数１１の関係を通し、エネルギー伝達量に関しても高い近似精度が得られると期待される。
【００６４】
ところで従来の技術では、ＳＧＳ応力のモデルとしてスマゴリンスキーモデル，ミックスドモデル，ツーパラメータ・ミックスドモデルが適用されている。これらのモデルはすべて、モデル中に渦粘性型のモデルを含んだものである。渦粘性型のモデルを含んだモデルを採用していることには理由がある。ミックスドモデルやツーパラメータ・ミックスドモデルの非渦粘性モデル部はスケール相似則モデルになっている。一般にスケール相似則モデルは、グリッドスケールからサブグリッドスケールへのエネルギー伝達量を過小評価する傾向があり、それ単独での解析では数値的に不安定になりやすい。そのため渦粘性係数を適切に選べば十分なエネルギー伝達量を与えられる渦粘性型のモデルと併用し、数値的に安定化させている。
【００６５】
渦粘性型のＳＧＳモデルでは、数２に示すように、ＳＧＳ応力の非等方成分テンソルが、グリッドスケールの歪率テンソルに比例すると仮定している。しかしこの仮定が成立しない場合も多いことが知られており、それは例えば東京大学出版会 数値流体力学シリーズ ３乱流解析 （１９９５年） 第８３頁などに記載されている。
【００６６】
つまり数値的な安定性を求めるためには渦粘性型のモデルが必須であるが、しかし渦粘性型のモデルではＳＧＳ応力テンソルの成分については十分には近似できないということになる。このような渦粘性型のモデルを含んだＳＧＳモデルに対し、リリーの最小自乗法を適用する場合、ＳＧＳ応力テンソルの成分に関して十分に高い近似精度を得ることは難しい。その結果、リリーの最小自乗法自体がエネルギー伝達特性の近似手法でないこともあり、グリッドスケールからサブグリッドスケールへのエネルギー伝達量に関しても、高い近似精度は期待できない。こうしたことにより従来の技術では、計算が数値的に不安定になりやすかったものと考えられる。
【００６７】
以上の考察より、上記の目的を達成するためには波数空間上でのエネルギー伝達特性を十分に反映させる必要があると考えられる。
【００６８】
ここで次の仮定を設ける。すなわち、
「テストスケールからサブテストスケールへのエネルギー伝達量は、グリッドスケールからサブグリッドスケールへのエネルギー伝達量と、サブテストスケールでありかつグリッドスケールである波数領域からのエネルギーの粘性散逸量との和に等しい」
と仮定する。図２を用いてこの仮定を説明する。図２は波数空間上でのエネルギースペクトルの概念図を示している。上記の仮定は、サブテストスケールでありかつグリッドスケールである波数領域７でのエネルギー授受の平衡を仮定するものである。つまりこの領域７に入ってくるエネルギーは、テストスケールからサブテストスケールへのエネルギー伝達量６であり、この領域７から出ていくエネルギーは、グリッドスケールからサブグリッドスケールへのエネルギー伝達量３′および、運動エネルギーのＳＴＳ成分でありかつＧＳ成分である成分の粘性散逸量５であって、こうしたエネルギー流入流出がバランスすると仮定する。この仮定は厳密には正しくないが、領域７の範囲がそれほど広くない場合には、近似的に成立すると考えられる。ＬＥＳの計算においては、領域７の範囲は、グリッドフィルタとテストフィルタとの二重フィルタの特性幅に対応する波数、およびグリッドフィルタの特性幅に対応する波数の間の範囲となる。これらの波数の比はおよそ２程度の値となるが、乱流場では乱流領域にあたる慣性小領域の範囲が数桁のオーダーに達することもあることを考えると、上記の仮定は決して悪くはない近似だと考えられる。
【００６９】
ちなみに運動エネルギーのＳＴＳ成分でありかつＧＳ成分である成分の粘性散逸量５は、運動エネルギーのＧＳ成分の粘性散逸量１′と、ＴＳ成分の粘性散逸量４との差で与えられる。
【００７０】
以降、上記の仮定を「局所スケール間平衡仮定」とよぶことにする。局所スケール間平衡仮定を表す式は、数１７で与えられる。
【００７１】
【数１７】
【００７２】
以降、数１７を「局所スケール間平衡仮定関係式」とよぶ。ここで運動エネルギーのＴＳ成分の粘性散逸量は数１８で与えられる。
【００７３】
【数１８】
【００７４】
またテストスケールからサブテストスケールへのエネルギー伝達量は数１９で与えられる。
【００７５】
【数１９】
【００７６】
局所スケール間平衡仮定関係式は、波数空間上でのエネルギー伝達特性をモデル化したものであり、グリッドスケールからサブグリッドスケールへのエネルギー伝達量自体をその式中に含むものである。従って局所スケール間平衡仮定関係式を満たすようにＳＧＳ応力モデルのモデル係数を調整してやれば、エネルギー伝達特性を十分に反映した、数値的に安定なＳＧＳモデル定数の自動算出法が構築できると考えられる。
【００７７】
本発明の実施の形態では、例えば、運動エネルギーのＴＳ成分の粘性散逸量の近似値と、運動エネルギーのＧＳ成分の粘性散逸量の近似値と、テストスケールからサブテストスケールへのエネルギー伝達量の近似値とを局所スケール間平衡仮定関係式に代入する。さらにグリッドスケールからサブグリッドスケールへのエネルギー伝達量を局所スケール間平衡仮定関係式に代入するに際し、グリッドスケールからサブグリッドスケールへのエネルギー伝達量を、ＳＧＳ応力テンソルと速度の空間変動のＧＳ成分に対する歪率テンソルとの内積に(−１)を乗じたのち、テストフィルタを作用させて得られる値とみなす。さらにこのときＳＧＳ応力テンソルを、その値をまだ定めていないモデル係数を１つだけ含むモデル式で与える。しかるのち、ＳＧＳ応力テンソルのモデル式中の、値がまだ定められていないモデル係数について局所スケール間平衡仮定関係式を整理することにより、未定であったモデル係数の値を算出するためのモデル係数算出式が求まるから、これを用いてモデル係数を算出する。以上の手続きを含んだＬＥＳの解析方法を提供することにより、上記の目的は達成される。
【００７８】
また、局所スケール間平衡仮定関係式のグリッドスケールからサブグリッドスケールへのエネルギー伝達量の評価に、その定義式である数１１の関係を用いることが望ましい。このときＳＧＳ応力テンソルのモデル式には、未定のモデル係数を１つだけ含むものを用いるが、これは局所スケール間平衡仮定関係式がスカラー式であるからで、未知量は１つのみでなければならないからである。
【００７９】
また、運動エネルギーのＴＳ成分およびＧＳ成分の粘性散逸量の差を０と近似して与えることも可能であり、対象とする乱流場の波数空間上での乱流領域のレンジが極めて広い場合など、粘性散逸の効果を無視できる場合に演算量を減らせる利点をもつ。
【００８０】
また、テストスケールからサブテストスケールへのエネルギー伝達量を、STS応力テンソルと速度の空間変動のＴＳ成分に対する歪率テンソルとの内積に、（−１）を乗じた値で与え、さらにこのとき、ＳＴＳ応力テンソルを、その値がまだ定められていないモデル係数を、１つも含まないモデル式で与えても良い。
【００８１】
テストスケールからサブテストスケールへのエネルギー伝達量を、その定義式である数１９で与えると、更にその解析精度が高まることが期待される。また、ＳＴＳ応力テンソルのモデルに、未定のモデル係数が１つも含まれていないから、局所スケール間平衡仮定関係式に含まれる未知数はＳＧＳ応力のモデル係数１つのみで、容易にモデル係数算出式を導出できる。
【００８２】
また、運動エネルギーのＧＳ成分の粘性散逸量を、速度の空間変動のＧＳ成分に対する歪率テンソルの絶対値の２乗と、動粘性係数との積の２倍の値に、さらにテストフィルタを作用させて得られた値、またはそれと等価の値で与える。また、運動エネルギーのＴＳ成分の粘性散逸量を、速度の空間変動のＴＳ成分に対する歪率テンソルの絶対値の２乗と、動粘性係数との積の２倍の値、もしくはそれと等価の値で与えても良い。運動エネルギーのＴＳ成分，ＧＳ成分の粘性散逸量を、その定義式である数１０，数１８で与えるので、更なる解析精度の向上が望める。
【００８３】
テストスケールからサブテストスケールへのエネルギー伝達量，運動エネルギーのＴＳ成分およびＧＳ成分の粘性散逸量に関し、物理的に正しい定義式を用いることにより、高い解析精度が得られるものと期待される。
【００８４】
ＳＧＳ応力を、渦粘性型のモデル式を含んだモデル式で与え、ＳＧＳ応力テンソルの近似において、スマゴリンスキーモデルやミックスドモデルなど、渦粘性型のモデル部を含むモデル式を適用しても良い。渦粘性型のモデルは、グリッドスケールからサブグリッドスケールへのエネルギー伝達量を近似するモデルとして適しており、計算の安定性向上にも寄与する。よって、解析精度の向上と、計算安定性の向上とが期待できる。
【００８５】
（実施例）
以下、本発明の実施例を図３により説明する。図３は、乱流の解析方法に、ＬＥＳを用いた場合の流れの非定常計算のフローチャートの一例を示すものである。また本例は密度一定の非圧縮性流れの場合を想定したものである。
【００８６】
まず初めに解析領域の形状データ，物性値，時間刻みといった、計算条件等の入力をＳ１において行う。ついでＳ２において速度，圧力等の初期値を設定する。ついでＳ３において、解析時刻ｔをゼロクリアする。ついでＳ４において、数２０に示す運動量保存式の、移流項による速度の加速量を計算し、計算機メモリに保存しておく。
【００８７】
【数２０】
【００８８】
ついでＳ５において、ＳＧＳ応力項および粘性項による速度の加速量を計算し、計算機メモリに保存しておく。ついでＳ６において、圧力項による速度の加速量を計算し、計算機メモリに保存しておく。ついでＳ７において、計算機メモリ上に保存された移流項，ＳＧＳ応力項および粘性項，圧力項による速度の加速量を用いて速度場を更新する。ついでＳ８においてＳ７で更新された速度場が、質量保存式を満足するか判定する。質量保存を十分には満足していない場合にはＳ９で速度・圧力の修正計算を行った後、再度Ｓ８で質量保存則の満足度をチェックする。Ｓ８において速度場が質量保存則を満たすと判断された場合は、Ｓ10において計算結果をモニタあるいは磁気ディスク上のファイルなどに出力する。ついでＳ１１において解析時刻ｔが、解析終了時刻ｔｅｎｄを越えたかどうか判定する。解析終了時刻を越えていない場合は、Ｓ１２において解析時刻を更新し、さらにＳ４にもどり、Ｓ４からＳ１０までの手続きを行ったのち、再度Ｓ１１の判定を行う。Ｓ１１で解析時刻が解析終了時刻を越えたと判断されたときは、計算を終了する。
【００８９】
本発明は図３のフローチャートでのＳ５のＳＧＳ応力項および粘性項の計算において、ＳＧＳ応力テンソルを算出するために用いる。
【００９０】
ここで本発明に係る第１の実施例につき説明する。
【００９１】
まずここで、グリッドフィルタをトップハットフィルタとみなし、さらにグリッドフィルタの特性幅について、その４乗以上に関する項を無視して、数２１で近似する。
【００９２】
【数２１】
【００９３】
つぎにテストフィルタを数２２で定義する。
【００９４】
【数２２】
【００９５】
このときテストフィルタの特性幅として、グリッドフィルタの特性幅の１乃至２倍程度の値を与えておく。
【００９６】
数２１および数２２を用いて数１３で定義されるＳＴＳ応力を整理し、グリッドフィルタの特性幅についてその４乗以上に関する項を無視すると、ＳＴＳ応力を近似するモデル式が数２３で与えられる。
【００９７】
【数２３】
【００９８】
数２３のモデル式は、値が未定であるモデル係数を、１つも含んでいない。
【００９９】
数１０，数１１，数１８，数１９，数２３を局所スケール間平衡仮定関係式である数１７に代入すると、数２４が得られる。
【０１００】
【数２４】
【０１０１】
ついでＳＧＳ応力テンソルを数２５のモデル式で近似する。数２５のＳＧＳ応力テンソルを近似するモデル式は、渦粘性モデル部と、非渦粘性モデル部からなる。非渦粘性モデル部は、その値が未定であるモデル係数を含まないものとする。このとき、数２５のモデル式に含まれる、値が未定のモデル係数は渦動粘性係数のみである。
【０１０２】
【数２５】
【０１０３】
数２５を数２４に代入して整理すると、数２６が得られる。
【０１０４】
【数２６】
【０１０５】
数２６は数２５の唯一のモデル係数である渦動粘性係数の算出式である。数２６の右辺は既知量のみから構成されるので、数２６を用いれば、渦動粘性係数を容易に求めることができる。
【０１０６】
数２６より算出した渦動粘性係数をもってＳＧＳ応力を近似し、図３のＳ５においてＳＧＳ応力項を計算すれば、流動状況に応じた、瞬時・局所的なモデル係数の最適値をもちいて解析したことになるので、高い解析精度が期待できる。さらに本実施例によれば、グリッドスケールからサブグリッドスケールへのエネルギー伝達特性が十分に反映されるため、数値的に安定であることが期待できる。
【０１０７】
ここで本発明に係る第２の実施例について説明する。ここで数２５のＳＧＳ応力モデルの非渦粘性モデル部を数２７で与える。
【０１０８】
【数２７】
【０１０９】
すなわちＳＧＳ応力を、渦粘性モデルのみで近似する。数２５，数２７を数２４に代入して整理すると、渦動粘性係数の算出式である数２８が得られる。
【０１１０】
【数２８】
【０１１１】
数２８より算出した渦動粘性係数をもってＳＧＳ応力を近似し、図３のＳ５においてＳＧＳ応力項を計算すれば、第１の実施例で述べた効果に加えて、ＳＧＳ応力のモデルの、非渦粘性モデル部については計算しなくてすむため、計算時間の短縮と使用メモリ量の節約が可能となる。
【０１１２】
ここで、本発明に係る第３の実施例について説明する。ここで数２５のＳＧＳ応力モデルの非渦粘性モデル部を、数２９で与える。
【０１１３】
【数２９】
【０１１４】
すなわち、非渦粘性モデル部を、バーディナモデルで与える。数２５，数２９を数２４に代入して整理すると、渦動粘性係数の算出式である数３０が得られる。
【０１１５】
【数３０】
【０１１６】
数３０より算出した渦動粘性係数をもってＳＧＳ応力を近似し、図３のＳ５においてＳＧＳ応力項を計算すれば、第１の実施例で述べた効果に加えて、バーディナモデル部がＳＧＳ応力テンソル成分について高い近似精度を有することからして、さらに高い精度の解析結果が得られると期待できる。
【０１１７】
ここで、本発明に係る第４の実施例について説明する。ここで数２５のＳＧＳ応力モデルの非渦粘性モデル部を、数３１で与える。
【０１１８】
【数３１】
【０１１９】
このとき、数３１に含まれるモデル定数の値は、あらかじめ設定しておく。数３１のモデル定数の設定については、まずＳＧＳ応力を数３２に示すモデルで与え、これに対して従来の技術で用いられてきたリリーの最小自乗法を適用し、算出する。
【０１２０】
【数３２】
【０１２１】
したがって数３１を数２５に代入して得られるＳＧＳ応力のモデル式では、値が未定のモデル係数はやはり、渦動粘性係数のみである。渦動粘性係数の算出式は、数２５，数３１を数２４に代入して整理することにより得られ、数３３で表される。
【０１２２】
【数３３】
【０１２３】
本実施例においても、第３の実施例と同様の効果が期待される。
【０１２４】
ここで、本発明に係る第５の実施例について説明する。ここで数２５のＳＧＳ応力モデルの非渦粘性モデル部を、数３４で与える。
【０１２５】
【数３４】
【０１２６】
数３４は、ＳＧＳ応力をテイラー展開を用いて近似したものである。数２５，数３４を数２４に代入して整理することにより、本実施例での唯一のモデル係数である、渦動粘性係数の算出式が数３５のように得られる。
【０１２７】
【数３５】
【０１２８】
本実施例においても、第３の実施例と同様の効果が期待される。
【０１２９】
ここで本発明に係る第６の実施例について説明する。ここで数２５のＳＧＳ応力モデルの渦動粘性係数を数４で与える。すなわち数２５のＳＧＳ応力モデルの渦粘性モデル部を、スマゴリンスキーモデルで与える。このとき、値が未定のモデル係数は、モデル定数Ｃのみである。数４，数２５を数２４に代入して整理すると、モデル定数Ｃの算出式である数３６が得られる。
【０１３０】
【数３６】
【０１３１】
数３６より算出したモデル定数ＣをもってＳＧＳ応力を近似し、図３のＳ５においてＳＧＳ応力項を計算すれば、第１の実施例で述べたのと同様の効果が期待できる。また、本実施例において、数２５で与えられるＳＧＳ応力モデルの非渦粘性モデル部については、第２，第３，第４，第５の実施例で述べたのと同様のものが適用できる。
【０１３２】
ここで本発明に係る第７の実施例について説明する。ここではまず、グリッドフィルタをトップハットフィルタとみなし、グリッドフィルタの特性幅の４乗に関する項まで考慮して、数３７で近似する。
【０１３３】
【数３７】
【０１３４】
つぎにテストフィルタを数２２で与え、２次中心差分法を用いて近似する。このとき、数２２の中の２階微分項は、打切り誤差について４次の項まで考慮すると、数３８のように表せる。
【０１３５】
【数３８】
【０１３６】
すなわち、数２２を２次中心差分法で離散化すると、グリッドフィルタの特性幅の４乗に関する打ち切り誤差項があらわれるため、テストフィルタとしては、実際には数３９のフィルタを用いたことになる。
【０１３７】
【数３９】
【０１３８】
数３７および数３９を用いて数１３で定義されるＳＴＳ応力を整理し、グリッドフィルタの特性幅についてその４乗以下に関する項までを考慮すると、ＳＴＳ応力を近似するモデル式が数４０で与えられる。
【０１３９】
【数４０】
【０１４０】
数４０のモデル式は、値が未定であるモデル係数を、１つも含んでいない。
【０１４１】
ＳＧＳ応力のモデル式である数２５と、ＳＴＳ応力のモデル式である数４０を局所スケール間平衡仮定関係式である数１７に代入して整理すると、数４１が得られる。
【０１４２】
【数４１】
【０１４３】
数４１は数２５の唯一のモデル係数である渦動粘性係数の算出式である。数４１の右辺は既知量のみから構成されるので、数４１を用いれば、渦動粘性係数を容易に求めることができる。
【０１４４】
数４１より算出した渦動粘性係数をもってＳＧＳ応力を近似し、図３のＳ５においてＳＧＳ応力項を計算すれば、流動状況に応じた、瞬時・局所的なモデル係数の最適値をもちいて解析したことになるので、高い解析精度が期待できる。さらに本実施例によれば、グリッドスケールからサブグリッドスケールへのエネルギー伝達特性が十分に反映されるため、数値的に安定であることが期待できる。さらに本実施例は、第１，第２，第３，第４，第５，第６の実施例がＳＴＳ応力をグリッドフィルタの特性幅の２乗に関する項まで考慮したのに対し、さらに加えて４乗に関する項につていも考慮しているため、ＳＴＳ応力に関して近似精度が高い。それにともないテストスケールからサブテストスケールへのエネルギー伝達量に関する近似精度も高くなるので、第１から第６の実施例に比べてさらに高精度の解析結果が期待できる。
【０１４５】
次に、各々のデータフローを説明する。図４は、ＳＧＳ応力テンソルの渦粘性モデル部における、渦動粘性係数を算出する際のデータフローを示すものである。ＧＳ成分の速度ベクトルデータ１１は、すでに求まっているものとする。１２は、ＧＳ成分の歪率テンソルの算出部であり、ＧＳ成分の速度ベクトル１１を用いてＧＳ成分の歪率テンソルを算出する。１３では、運動エネルギーのＧＳ成分の粘性散逸量を、１２で求めたＧＳ成分の歪率テンソルを用いて算出する。１４は、ＳＧＳ応力テンソルの非渦粘性モデル部の算出部であり、ＧＳ成分の速度ベクトル１１を用いて、例えば数２７，数２９，数３１，数３４等の近似式により算出する。１５では、ＴＳ成分の速度ベクトルをＧＳ成分の速度ベクトル１１を用いて算出する。算出に用いるテストフィルタには、例えば数２２，数３９等を適用する。１６はＴＳ成分の歪率テンソルの算出部であり、１５で求めたＴＳ成分の速度ベクトルを用いて算出する。１７では、運動エネルギーのＴＳ成分の粘性散逸量を、１６で求めたＴＳ成分の歪率テンソルを用いて算出する。１８はＳＴＳ応力テンソルの算出部で、１１のＧＳ成分の速度ベクトル、１５のＴＳ成分の速度ベクトル等を用い、例えば数２３，数４０等の近似式によってＳＴＳ応力テンソルの近似値を算出する。１９では渦動粘性係数を算出する。
【０１４６】
算出に際しては、１２のＧＳ成分の歪率テンソル、１６のＴＳ成分の歪率テンソル、１３の運動エネルギーのＧＳ成分の粘性散逸量、１７の運動エネルギーのＴＳ成分の粘性散逸量、１４のＳＧＳ応力の非渦粘性モデル部、１８のＳＴＳ応力テンソルを用いる。なお、ここで求めた渦動粘性係数は、ＳＧＳ応力テンソルの近似式である数２あるいは数２５を介し、フローチャート図３のステップＳ５において、ＳＧＳ応力項の計算で用いられる。
【０１４７】
また、図５はＳＧＳ応力テンソルの渦粘性モデル部における、渦動粘性係数を算出する際のデータフローを示すものであり、ＴＳ成分の歪率テンソル算出部１６′以外は、図４と同様である。１６′はＴＳ成分の歪率テンソルの算出部であり、１２のＧＳ成分の歪率テンソルに対し、テストフィルタ操作を行うことにより算出する（図４ではＴＳ成分の速度ベクトル１５を用いて算出）。
【０１４８】
また、図６は、ＳＧＳ応力テンソルの渦粘性モデル部における、渦動粘性係数に関するモデル定数を算出する際のデータフローを示すものである。基本的には図４のデータフローと同様である。両者の違いは、図４が渦動粘性係数を求めるデータフローであったのに対し、図６が渦粘性モデルとしてスマゴリンスキーモデルを用いる際（数４）の、モデル定数Ｃを算出するデータフローである点である。２０では渦動粘性係数に関するモデル定数を算出する。算出に際しては、１２のＧＳ成分の歪率テンソル、１６のＴＳ成分の歪率テンソル、１３の運動エネルギーのＧＳ成分の粘性散逸量、１７の運動エネルギーのＴＳ成分の粘性散逸量、１４のＳＧＳ応力の非渦粘性モデル部、１８のＳＴＳ応力テンソルを用いる。
【０１４９】
また、図７は、ＴＳ成分の歪率テンソル算出部１６′以外は、図６と同様である。
【０１５０】
以上のように、本実施例では、乱流の解析方法の一つであるＬＥＳにおいて、ＳＧＳ応力のモデル式に含まれるモデル係数の局所的な値を、瞬時の流動状況に応じて自動的に算出するため、モデル係数のチューニングを必要とせずに、高い精度の解析結果が得られるという効果を奏する。
【０１５１】
ここで、この乱流解析方法によるプログラムを作成し、このプログラム若しくはこのプログラムを記録した記録媒体でコンピュータによって容易に乱流解析を実施できる。
【０１５２】
またＬＥＳのＳＧＳ応力のモデル係数を算出するに際し、本実施例では波数空間上でのエネルギー伝達特性のモデル式である、局所スケール間平衡仮定関係式をもちいる。さらに本実施例では特に、グリッドスケールからサブグリッドスケールへのエネルギー伝達量の評価に、その定義式である数１１の関係を用いる。そのため本実施例によれば、ＳＧＳ応力のモデル係数の算出に際し、エネルギー伝達特性が十分に考慮されるから、解析が数値的に安定化するという効果が得られる。
【０１５３】
更に、テストスケールからサブテストスケールへのエネルギー伝達量を、その定義式である数１９で与えるので、高い解析精度が得られるという効果がある。
【０１５４】
更に、運動エネルギーのＴＳ成分，ＧＳ成分の粘性散逸量を、その定義式である数１０，数１８で与えるので、さらに高い解析制度が得られるという効果がある。
【０１５５】
更に、テストスケールからサブテストスケールへのエネルギー伝達量，運動エネルギーのＴＳ成分およびＧＳ成分の粘性散逸量に関し、物理的に正しい定義式を用いたことにより、高い解析精度が得られる効果がある。
【０１５６】
更に、ＳＧＳ応力を、渦粘性型のモデル式を含んだモデル式で与えるものである。渦粘性型のモデルは計算の安定化に寄与するから、本実施例によれば、数値的に極めて安定な乱流解析が可能となるという効果が得られる。
【０１５７】
【発明の効果】
本発明によると、乱流の解析方法の一つであるＬＥＳにおいて、ＳＧＳ応力のモデル式中に含まれるモデル係数もしくはモデル定数を、自動的に算出し、数値的に安定な乱流の解析方法，乱流の解析システム並びに乱流解析プログラムを提供することができるという効果を奏する。
【図面の簡単な説明】
【図１】波数空間上での運動エネルギー伝達の概念図。
【図２】局所スケール間平衡仮定の説明図。
【図３】ＬＥＳをもちいた流れの非定常計算のフローチャート。
【図４】渦動粘性係数の算出に関するデータフロー。
【図５】渦動粘性係数の算出に関するデータフロー。
【図６】渦動粘性係数に関するモデル定数の算出に関するデータフロー。
【図７】渦動粘性係数に関するモデル定数の算出に関するデータフロー。
【符号の説明】
１…運動エネルギーのＧＳ成分の粘性散逸量、２…運動エネルギーのＳＧＳ成分の粘性散逸量、３…グリッドスケールからサブグリッドスケールへのエネルギー伝達量、１′…運動エネルギーのＧＳ成分の粘性散逸量、２′…運動エネルギーのＳＧＳ成分の粘性散逸量、３′…グリッドスケールからサブグリッドスケールへのエネルギー伝達量、４…運動エネルギーのＴＳ成分の粘性散逸量、５…運動エネルギーのＳＴＳ成分でありかつＧＳ成分である成分の粘性散逸量、６…テストスケールからサブテストスケールへのエネルギー伝達量、７…サブテストスケールでありかつグリッドスケールである波数領域、１１…ＧＳ成分の速度ベクトルのデータ、１２…ＧＳ成分の歪率テンソルの算出部、１３…運動エネルギーのＧＳ成分の粘性散逸量の算出部、１４…ＳＧＳ応力テンソルの非渦粘性モデル部の算出部、１５…ＴＳ成分の速度ベクトルの算出部、１６…ＴＳ成分の歪率テンソルの算出部、１６′…ＴＳ成分の歪率テンソルの算出部、１７…運動エネルギーのＴＳ成分の粘性散逸量の算出部、１８…ＳＴＳ応力テンソルの算出部、１９…渦動粘性係数の算出部、２０…渦動粘性係数に関するモデル定数の算出部、Ｓ１…計算条件等の入力ステップ、Ｓ２…初期値設定ステップ、Ｓ３…解析時刻ゼロクリアステップ、Ｓ４…移流項の計算ステップ、Ｓ５…ＳＧＳ応力項・粘性項の計算ステップ、Ｓ６…圧力項の計算ステップ、Ｓ７…速度場の更新ステップ、Ｓ８…質量保存則の満足度の判定ステップ、Ｓ９…速度・圧力の修正計算ステップ、Ｓ１０…計算結果の出力ステップ、Ｓ１１…計算終了の判定ステップ、Ｓ１２…解析時刻の更新ステップ。
【発明の属する技術分野】
本発明は流体機器や流体関連機器における流れの数値解析方法に関わり、特に乱流場を非定常解析する計算方法であるラージ・エディ・シミュレーションにおいて、ＳＧＳ応力のモデル係数を自動的に決定する方法に関する。
【０００２】
【従来の技術】
流体機器の設計や開発では、流れの数値解析が利用されることが多い。そして対象とする流れ場の多くは乱流である。乱流の解析方法の一つにラージ・エディ・シミュレーション（ＬＥＳ）がある。
【０００３】
流れの数値解析手法の多くは、解析の対象となる領域を、空間的に有限個の領域、あるいは有限個の物理量定義点に分割して扱う。これを空間の離散化といい、一般的にはメッシュ分割，要素分割，セル分割などとよぶ。空間を離散化すると、たとえばメッシュ幅よりもサイズが小さい渦などの流れ構造は、もはや解像できなくなってしまう。しかし、こうした微小な流れ構造も、全体の流れ場に影響を及ぼすので、その影響を考慮する必要がある。こうした微小な流れ構造（波数空間上での高波数域の速度の変動成分）の影響をサブグリッドスケールモデル（ＳＧＳモデル）として考慮し、流れを非定常計算するのがＬＥＳである。
【０００４】
空間の離散化は、物理量の空間的な変動の、高波数成分を除去する一種のローパスフィルタとみなすことができ、このフィルタをグリッドフィルタとよぶ。波数空間において、グリッドフィルタにより除去される高波数の領域のことをサブグリッドスケールとよび、除去されない低波数域のことをグリッドスケールとよぶ。また波数空間上でグリッドスケールに含まれる物理量の空間変動成分のことをグリッドスケール成分（ＧＳ成分）とよび、サブグリッドスケールに含まれる空間変動成分をサブグリッドスケール成分（ＳＧＳ成分）とよぶ。
【０００５】
任意の物理量をφであらわすとき、空間的に分布するφにグリッドフィルタを作用させて得られる値を
【０００６】
【数４２】
であらわす。
【０００８】
空間の離散化で解像できなくなった微小な流れ構造、すなわち空間的な速度変動のＳＧＳ成分が、流れの運動量変化に及ぼす影響は、応力の形であらわれる。この応力をサブグリッドスケール応力（ＳＧＳ応力）という。
【０００９】
速度のｉ方向成分をｕi と表すとき、ＳＧＳ応力は数１であらわされる。
【００１０】
【数１】
ＳＧＳ応力中の
【００１２】
【数４３】
については直接的には計算することができないため、ＳＧＳ応力に対してモデル化が必要となる。
【００１４】
ＳＧＳ応力に対する代表的なＳＧＳモデルに、スマゴリンスキーモデルがある。スマゴリンスキーモデルは数２で表される渦粘性型のＳＧＳモデルである。
【００１５】
【数２】
ここでグリッドスケールでの歪率テンソルは数３で与えられる。
【００１７】
【数３】
スマゴリンスキーモデルでは渦動粘性係数をモデル定数Ｃを含んだ数４で与える。
【００１９】
【数４】
なお、本明細書においては数式の記述において、アインシュタインの総和規約を用いるものとする。また任意のテンソルΦijに対し、Φij * はその非等方成分Φij *≡Φij−δijΦkk／３をあらわすものとする。
【００２１】
スマゴリンスキーモデルのモデル定数は流れの性状によりその最適値が異なる。例えば一様等方性乱流ではその最適値は０.０３ から０.０４ 程度の値であり、チャンネル乱流では０.０１ 程度の値が最適値である。
【００２２】
そのためスマゴリンスキーモデルを用いた解析で精度の高い結果を得るためには、対象とする流れ場に適したモデル定数を用いる必要がある。しかし複雑な流れ場を対象とするとき、そのようなモデル定数のチューニングは極めて難しい。
【００２３】
そこで流動状態に応じてモデル定数の局所的な値を自動的に算出する方法が提案されている。その方法はフィジックス オブ フルイズ エー ４（１９９２年）第６３３頁から第６３５頁（Physics of Fluids A 4 (1992) pp.633-635 ）に記載されたものである。この方法をもちいたＳＧＳモデルはダイナミック・スマゴリンスキーモデル（ＤＳＭ）と呼ばれる。またこのときモデル定数を算出するために用いる計算手続きはリリーの最小自乗法とよばれる。
【００２４】
ＳＧＳモデルにスマゴリンスキーモデル以外のモデルを適用する場合についても、モデル定数を自動的に算出する方法が提案されている。例えば数５はスケール相似則モデルの一つであるバーディナモデルであるが、これと数２，数４のスマゴリンスキーモデルとを線形結合したモデルとして、数６に示すミックスドモデルがある。
【００２５】
【数５】
【数６】
数６のミックスドモデルを対象として、式中に含まれるスマゴリンスキーモデル部のモデル定数Ｃを算出する方法が、フィジックス オブ フルイズ Ａ ５（１９９３年）第３１８６頁から第３１９６頁（Physics of Fluids A 5(1993)pp.3186-3196）に記載されている。この方法はダイナミック・ミックスドモデル（ＤＭＭ）とよばれる。ダイナミック・ミックスドモデルにおいても、モデル定数の算出手続きにはリリーの最小自乗法が用いられている。
【００２８】
また数７に示すツーパラメータ・ミックスドモデルはミックスドモデルの一種であり、モデル定数を２つ含むモデルである。
【００２９】
【数７】
この２つのモデル定数を自動的に算出する方法はフィジックス オブ フルイズ ７（１９９５年）第２８３１頁から第２８４７頁（Physics of Fluids 7 (1997) pp.2831-2847)に記載されている。この方法はダイナミック・ツーパラメータ・ミックスドモデル（ＤＴＭＭ）と呼ばれる。またこの方法でモデル定数を算出するために用いられる方法はやはり、リリーの最小自乗法である。
【００３１】
【発明が解決しようとする課題】
上記の従来技術における問題点は、数値的に不安定な点である。すなわち、ＳＧＳモデルのモデル定数の自動算出に従来の技術を用いた解析では、数値振動が生じやすく、場合によっては非物理的な振動解が計算結果として得られる。さらには計算途中においてゼロ割りやフローティングポイントオーバーフローなどが生じ、計算の続行が不能となることもある。
【００３２】
こうした数値不安定を防ぐ目的で、統計的に一様な方向への平均操作を用いる方法が提案されており、その一例はフィジックス オブ フルイズ エー ３（１９９１年）第１７６０頁から第１７６５頁（Physics of Fluids A 3 (1991)pp.1760-1765）に記載されている。但しこうした平均操作は単純な体系にしか適用できず、またその物理的な根拠も不明であるため、一般的な解決策ではない。乱流の解析方法の一つであるＬＥＳにおいて、ＳＧＳ応力のモデル式中に含まれるモデル係数もしくはモデル定数の局所的な値を、瞬時の流動状況に応じて自動的に算出することが望ましい。さらに加えて、数値的な不安定現象が生じにくい、数値的に安定な計算方法であることが望ましい。
【００３３】
本発明の目的は、乱流の解析方法の一つであるＬＥＳにおいて、ＳＧＳ応力のモデル式中に含まれるモデル係数もしくはモデル定数を、自動的に算出し、数値的に安定な乱流の解析方法，乱流の解析システム並びに乱流解析プログラムを提供することにある。
【００３４】
【課題を解決するための手段】
本発明は、該離散化に伴い空間的には解像することができなくなった高波数域の速度の空間変動成分が流体の運動に与える影響であるサブグリッドスケール応力（ＳＧＳ応力）をモデル化したモデル式を計算する計算手順が前記計算機メモリもしくは記録媒体に保存され、物理量の空間変動の波数成分のうち該離散化に伴い空間的には解像することができなくなった高波数域 ( サブグリッドスケール ) の成分をサブグリッドスケール成分（ＳＧＳ成分），物理量の空間変動の波数成分のうちＳＧＳ成分ではない低波数域（グリッドスケール）の成分をグリッドスケール成分（ＧＳ成分）とし、物理量の空間変動の高波数成分を除去する作用をもつ空間フィルタであるテストフィルタを、物理量の空間変動の前記
ＧＳ成分に対して作用させ、前記ＧＳ成分のうちの高波数域の成分を除去した残りの低波数域（テストスケール）の成分をテストスケール成分（ＴＳ成分）とし、物理量の空間変動の波数成分のうち前記ＴＳ成分には含まれない高波数域（サブテストスケール）の成分をサブテストスケール成分（ＳＴＳ成分）とする場合、速度ベクトルの前記ＧＳ成分を算出し、前記計算メモリもしくは前記記録媒体に記憶する手段と、前記算出された前記速度ベクトルのＧＳ成分のデータに基づき、運動エネルギーのＴＳ成分の粘性散逸量の近似値と、運動エネルギーのＧＳ成分の粘性散逸量の近似値と、テストスケールからサブテストスケールへのエネルギー伝達量の近似値とを算出する手段と、運動エネルギーの該ＴＳ成分の粘性散逸量とテストスケールからサブテストスケールへのエネルギー伝達量との和が、運動エネルギーの該ＧＳ成分の粘性散逸量とグリッドスケールからサブグリッドスケールへのエネルギー伝達量との和に等しくなるという、前記計算メモリもしくは前記記録媒体に保存された局所スケール間平衡仮定関係式の計算手順を呼び出す手段と、グリッドスケールからサブグリッドスケールへのエネルギー伝達量を、ＳＧＳ応力テンソルと、速度の空間変動のＧＳ成分に対する歪率テンソルとの内積に（−１）を乗じたのち、テストフィルタを作用させて得られる値であるとみなし、このとき該ＳＧＳ応力テンソルを、その値をまだ定めていないモデル係数を、１つだけ含む前記モデル式で与えるとき、該局所スケール間平衡仮定関係式を、該値をまだ定めていないモデル係数について整理することにより、該ＳＧＳ応力のモデル式に含まれる該モデル係数を算出するためのモデル係数算出式を導き、該モデル係数算出式を用いて該モデル係数を計算してグリッドスケールからサブグリッドスケールへのエネルギー伝達量を計算する手段とを備えたことを特徴とする。
【００３５】
【発明の実施の形態】
ＳＧＳ応力は流れの運動量ベクトルを時間・空間的に変化させる作用をもつ。そのため必然的に、流れの運動エネルギー変化に対しても作用を及ぼすことになる。
【００３６】
数８にサブグリッドスケール運動エネルギー（ＳＧＳ運動エネルギー）の時間変化に関する簡略式を示す。
【００３７】
【数８】
ＳＧＳ運動エネルギーは、速度の空間変動のＳＧＳ成分がもつ運動エネルギーに相当する。数８の中の保存項は、移流項や拡散項など、ＳＧＳ運動エネルギーの輸送に関連した項であり、ＳＧＳ運動エネルギーの生成や消滅には直接的には関係しない項である。ＳＧＳ運動エネルギーの生成・消滅に関与するのは、運動エネルギーのＳＧＳ成分の粘性散逸量と、グリッドスケールからサブグリッドスケールへのエネルギー伝達量である。数８において、粘性散逸は消滅項として働き、グリッドスケールからサブグリッドスケールへのエネルギー伝達は生成項として働く。運動エネルギーのＳＧＳ成分の粘性散逸量は数９で与えられ、運動エネルギーのＧＳ成分の粘性散逸量は数１０で与えられる。
【００３９】
【数９】
【数１０】
また、グリッドスケールからサブグリッドスケールへのエネルギー伝達量は数１１で与えれれる。
【００４２】
【数１１】
数１１に示されるように、ＳＧＳ応力が関与するのは、グリッドスケールからサブグリッドスケールへのエネルギー伝達量に関してである。
【００４４】
数１２にグリッドスケール運動エネルギー（ＧＳ運動エネルギー）の時間変化に関する簡略式を示す。
【００４５】
【数１２】
ＧＳ運動エネルギーは、速度の空間変動のＧＳ成分がもつ運動エネルギーに相当する。数１２において、運動エネルギーのＧＳ成分の粘性散逸は消滅項として働き、グリッドスケールからサブグリッドスケールへのエネルギー伝達量もまた消滅項として働く。すなわち、ＧＳ運動エネルギーの一部が消滅し、その消滅分がＳＧＳ運動エネルギーの生成分となる。
【００４７】
図１に、上に述べた波数空間での運動エネルギー伝達に関する概念図を示す。図１は波数空間上でのエネルギースペクトルを表しており、運動エネルギーの
ＧＳ成分の粘性散逸１、およびＳＧＳ成分の粘性散逸２により、運動エネルギーの一部は熱として散逸する。また、運動エネルギーのＧＳ成分の一部は、グリッドスケールからサブグリッドスケールへのエネルギー伝達３により、運動エネルギーのＳＧＳ成分に移動する。
【００４８】
上述した波数空間上での運動エネルギーの伝達特性は、計算の数値安定性および計算精度の観点から極めて重要である。グリッドスケールからサブグリッドスケールへのエネルギー伝達量を過小評価するということは、ＧＳ運動エネルギーの減衰が小さいということで、グリッドスケールの高波数域に運動エネルギーが蓄積されることになる。すなわち、非物理的な高波数の数値振動が生じやすくなる。その反対にエネルギー伝達量を過大評価すると、ＧＳ運動エネルギーの減衰が大きくなり、特に高波数域での減衰が大きくなる。そのため微小な流れ構造を解像することができなくなってしまい、解析結果の精度が悪くなる。したがってＳＧＳ応力のモデルには、グリッドスケールからサブグリッドスケールへのエネルギー伝達量を適切に模擬することが求められる。
【００４９】
従来の技術に記載したＳＧＳモデルのモデル定数の決定法はすべて、リリーの最小自乗法を用いたものであった。ここでリリーの最小自乗法について説明する。
【００５０】
まず、物理量の空間変動の高波数成分を除去する空間フィルタであるテストフィルタを導入する。このとき任意の物理量φにテストフィルタを作用させて得られる値を
【００５１】
【数４４】
と表記することにする。また、テストフィルタを物理量のＧＳ成分に作用させ、ＧＳ成分のうちの高波数域の成分を除去した、残りの低波数領域をテストスケールとよび、その成分をテストスケール成分（ＴＳ成分）とよぶ。また全波数域からテストスケールを除去した残りの高波数領域をサブテストスケールとよび、その成分をサブテストスケール成分（ＳＴＳ成分）とよぶことにする。
【００５３】
このとき、空間的な速度変動のＳＴＳ成分が運動量に及ぼす影響であるサブテストスケール応力（ＳＴＳ応力）は、数１３で与えられる。
【００５４】
【数１３】
ＳＧＳ応力とＳＴＳ応力の間には、数１４に示す関係式が成り立つ。
【００５６】
【数１４】
この関係式は恒等式であり、ジャーマノ恒等式と呼ばれる。
【００５８】
リリーの最小自乗法では、ジャーマノ恒等式中のＳＧＳ応力およびＳＴＳ応力に対し、同一のＳＧＳモデルを適用する。また、このときの残差テンソルを数１５で定義する。
【００５９】
【数１５】
ＳＧＳ応力に対するモデル式にＣ1からＣnまでのｎ個のモデル定数が含まれる場合について考える。リリーの最小自乗法では、数１６に示すように、残差テンソル成分の自乗和の、モデル定数に対する偏微分値が、すべてのモデル定数に対して０となるよう、モデル定数Ｃ1〜Ｃnを求める。
【００６１】
【数１６】
この手続きは一般の最小自乗法と同じである。すなわち、リリーの最小自乗法は、ジャーマノ恒等式をなるべく満足するように、モデル式で与えられたＳＧＳ応力テンソルの成分値を最適化する方法である。
【００６３】
このように、リリーの最小自乗法自体は、グリッドスケールからサブグリッドスケールへのエネルギー伝達特性に関しては、直接的には何一つ言及していない。しかし、ジャーマノ恒等式に関する残差テンソルが十分に小さくなるならば、ＳＧＳ応力テンソルの成分に関する近似精度が十分に高くなることが期待でき、数１１の関係を通し、エネルギー伝達量に関しても高い近似精度が得られると期待される。
【００６４】
ところで従来の技術では、ＳＧＳ応力のモデルとしてスマゴリンスキーモデル，ミックスドモデル，ツーパラメータ・ミックスドモデルが適用されている。これらのモデルはすべて、モデル中に渦粘性型のモデルを含んだものである。渦粘性型のモデルを含んだモデルを採用していることには理由がある。ミックスドモデルやツーパラメータ・ミックスドモデルの非渦粘性モデル部はスケール相似則モデルになっている。一般にスケール相似則モデルは、グリッドスケールからサブグリッドスケールへのエネルギー伝達量を過小評価する傾向があり、それ単独での解析では数値的に不安定になりやすい。そのため渦粘性係数を適切に選べば十分なエネルギー伝達量を与えられる渦粘性型のモデルと併用し、数値的に安定化させている。
【００６５】
渦粘性型のＳＧＳモデルでは、数２に示すように、ＳＧＳ応力の非等方成分テンソルが、グリッドスケールの歪率テンソルに比例すると仮定している。しかしこの仮定が成立しない場合も多いことが知られており、それは例えば東京大学出版会 数値流体力学シリーズ ３乱流解析 （１９９５年） 第８３頁などに記載されている。
【００６６】
つまり数値的な安定性を求めるためには渦粘性型のモデルが必須であるが、しかし渦粘性型のモデルではＳＧＳ応力テンソルの成分については十分には近似できないということになる。このような渦粘性型のモデルを含んだＳＧＳモデルに対し、リリーの最小自乗法を適用する場合、ＳＧＳ応力テンソルの成分に関して十分に高い近似精度を得ることは難しい。その結果、リリーの最小自乗法自体がエネルギー伝達特性の近似手法でないこともあり、グリッドスケールからサブグリッドスケールへのエネルギー伝達量に関しても、高い近似精度は期待できない。こうしたことにより従来の技術では、計算が数値的に不安定になりやすかったものと考えられる。
【００６７】
以上の考察より、上記の目的を達成するためには波数空間上でのエネルギー伝達特性を十分に反映させる必要があると考えられる。
【００６８】
ここで次の仮定を設ける。すなわち、
「テストスケールからサブテストスケールへのエネルギー伝達量は、グリッドスケールからサブグリッドスケールへのエネルギー伝達量と、サブテストスケールでありかつグリッドスケールである波数領域からのエネルギーの粘性散逸量との和に等しい」
と仮定する。図２を用いてこの仮定を説明する。図２は波数空間上でのエネルギースペクトルの概念図を示している。上記の仮定は、サブテストスケールでありかつグリッドスケールである波数領域７でのエネルギー授受の平衡を仮定するものである。つまりこの領域７に入ってくるエネルギーは、テストスケールからサブテストスケールへのエネルギー伝達量６であり、この領域７から出ていくエネルギーは、グリッドスケールからサブグリッドスケールへのエネルギー伝達量３′および、運動エネルギーのＳＴＳ成分でありかつＧＳ成分である成分の粘性散逸量５であって、こうしたエネルギー流入流出がバランスすると仮定する。この仮定は厳密には正しくないが、領域７の範囲がそれほど広くない場合には、近似的に成立すると考えられる。ＬＥＳの計算においては、領域７の範囲は、グリッドフィルタとテストフィルタとの二重フィルタの特性幅に対応する波数、およびグリッドフィルタの特性幅に対応する波数の間の範囲となる。これらの波数の比はおよそ２程度の値となるが、乱流場では乱流領域にあたる慣性小領域の範囲が数桁のオーダーに達することもあることを考えると、上記の仮定は決して悪くはない近似だと考えられる。
【００６９】
ちなみに運動エネルギーのＳＴＳ成分でありかつＧＳ成分である成分の粘性散逸量５は、運動エネルギーのＧＳ成分の粘性散逸量１′と、ＴＳ成分の粘性散逸量４との差で与えられる。
【００７０】
以降、上記の仮定を「局所スケール間平衡仮定」とよぶことにする。局所スケール間平衡仮定を表す式は、数１７で与えられる。
【００７１】
【数１７】
以降、数１７を「局所スケール間平衡仮定関係式」とよぶ。ここで運動エネルギーのＴＳ成分の粘性散逸量は数１８で与えられる。
【００７３】
【数１８】
またテストスケールからサブテストスケールへのエネルギー伝達量は数１９で与えられる。
【００７５】
【数１９】
局所スケール間平衡仮定関係式は、波数空間上でのエネルギー伝達特性をモデル化したものであり、グリッドスケールからサブグリッドスケールへのエネルギー伝達量自体をその式中に含むものである。従って局所スケール間平衡仮定関係式を満たすようにＳＧＳ応力モデルのモデル係数を調整してやれば、エネルギー伝達特性を十分に反映した、数値的に安定なＳＧＳモデル定数の自動算出法が構築できると考えられる。
【００７７】
本発明の実施の形態では、例えば、運動エネルギーのＴＳ成分の粘性散逸量の近似値と、運動エネルギーのＧＳ成分の粘性散逸量の近似値と、テストスケールからサブテストスケールへのエネルギー伝達量の近似値とを局所スケール間平衡仮定関係式に代入する。さらにグリッドスケールからサブグリッドスケールへのエネルギー伝達量を局所スケール間平衡仮定関係式に代入するに際し、グリッドスケールからサブグリッドスケールへのエネルギー伝達量を、ＳＧＳ応力テンソルと速度の空間変動のＧＳ成分に対する歪率テンソルとの内積に(−１)を乗じたのち、テストフィルタを作用させて得られる値とみなす。さらにこのときＳＧＳ応力テンソルを、その値をまだ定めていないモデル係数を１つだけ含むモデル式で与える。しかるのち、ＳＧＳ応力テンソルのモデル式中の、値がまだ定められていないモデル係数について局所スケール間平衡仮定関係式を整理することにより、未定であったモデル係数の値を算出するためのモデル係数算出式が求まるから、これを用いてモデル係数を算出する。以上の手続きを含んだＬＥＳの解析方法を提供することにより、上記の目的は達成される。
【００７８】
また、局所スケール間平衡仮定関係式のグリッドスケールからサブグリッドスケールへのエネルギー伝達量の評価に、その定義式である数１１の関係を用いることが望ましい。このときＳＧＳ応力テンソルのモデル式には、未定のモデル係数を１つだけ含むものを用いるが、これは局所スケール間平衡仮定関係式がスカラー式であるからで、未知量は１つのみでなければならないからである。
【００７９】
また、運動エネルギーのＴＳ成分およびＧＳ成分の粘性散逸量の差を０と近似して与えることも可能であり、対象とする乱流場の波数空間上での乱流領域のレンジが極めて広い場合など、粘性散逸の効果を無視できる場合に演算量を減らせる利点をもつ。
【００８０】
また、テストスケールからサブテストスケールへのエネルギー伝達量を、STS応力テンソルと速度の空間変動のＴＳ成分に対する歪率テンソルとの内積に、（−１）を乗じた値で与え、さらにこのとき、ＳＴＳ応力テンソルを、その値がまだ定められていないモデル係数を、１つも含まないモデル式で与えても良い。
【００８１】
テストスケールからサブテストスケールへのエネルギー伝達量を、その定義式である数１９で与えると、更にその解析精度が高まることが期待される。また、ＳＴＳ応力テンソルのモデルに、未定のモデル係数が１つも含まれていないから、局所スケール間平衡仮定関係式に含まれる未知数はＳＧＳ応力のモデル係数１つのみで、容易にモデル係数算出式を導出できる。
【００８２】
また、運動エネルギーのＧＳ成分の粘性散逸量を、速度の空間変動のＧＳ成分に対する歪率テンソルの絶対値の２乗と、動粘性係数との積の２倍の値に、さらにテストフィルタを作用させて得られた値、またはそれと等価の値で与える。また、運動エネルギーのＴＳ成分の粘性散逸量を、速度の空間変動のＴＳ成分に対する歪率テンソルの絶対値の２乗と、動粘性係数との積の２倍の値、もしくはそれと等価の値で与えても良い。運動エネルギーのＴＳ成分，ＧＳ成分の粘性散逸量を、その定義式である数１０，数１８で与えるので、更なる解析精度の向上が望める。
【００８３】
テストスケールからサブテストスケールへのエネルギー伝達量，運動エネルギーのＴＳ成分およびＧＳ成分の粘性散逸量に関し、物理的に正しい定義式を用いることにより、高い解析精度が得られるものと期待される。
【００８４】
ＳＧＳ応力を、渦粘性型のモデル式を含んだモデル式で与え、ＳＧＳ応力テンソルの近似において、スマゴリンスキーモデルやミックスドモデルなど、渦粘性型のモデル部を含むモデル式を適用しても良い。渦粘性型のモデルは、グリッドスケールからサブグリッドスケールへのエネルギー伝達量を近似するモデルとして適しており、計算の安定性向上にも寄与する。よって、解析精度の向上と、計算安定性の向上とが期待できる。
【００８５】
（実施例）
以下、本発明の実施例を図３により説明する。図３は、乱流の解析方法に、ＬＥＳを用いた場合の流れの非定常計算のフローチャートの一例を示すものである。また本例は密度一定の非圧縮性流れの場合を想定したものである。
【００８６】
まず初めに解析領域の形状データ，物性値，時間刻みといった、計算条件等の入力をＳ１において行う。ついでＳ２において速度，圧力等の初期値を設定する。ついでＳ３において、解析時刻ｔをゼロクリアする。ついでＳ４において、数２０に示す運動量保存式の、移流項による速度の加速量を計算し、計算機メモリに保存しておく。
【００８７】
【数２０】
ついでＳ５において、ＳＧＳ応力項および粘性項による速度の加速量を計算し、計算機メモリに保存しておく。ついでＳ６において、圧力項による速度の加速量を計算し、計算機メモリに保存しておく。ついでＳ７において、計算機メモリ上に保存された移流項，ＳＧＳ応力項および粘性項，圧力項による速度の加速量を用いて速度場を更新する。ついでＳ８においてＳ７で更新された速度場が、質量保存式を満足するか判定する。質量保存を十分には満足していない場合にはＳ９で速度・圧力の修正計算を行った後、再度Ｓ８で質量保存則の満足度をチェックする。Ｓ８において速度場が質量保存則を満たすと判断された場合は、Ｓ10において計算結果をモニタあるいは磁気ディスク上のファイルなどに出力する。ついでＳ１１において解析時刻ｔが、解析終了時刻ｔｅｎｄを越えたかどうか判定する。解析終了時刻を越えていない場合は、Ｓ１２において解析時刻を更新し、さらにＳ４にもどり、Ｓ４からＳ１０までの手続きを行ったのち、再度Ｓ１１の判定を行う。Ｓ１１で解析時刻が解析終了時刻を越えたと判断されたときは、計算を終了する。
【００８９】
本発明は図３のフローチャートでのＳ５のＳＧＳ応力項および粘性項の計算において、ＳＧＳ応力テンソルを算出するために用いる。
【００９０】
ここで本発明に係る第１の実施例につき説明する。
【００９１】
まずここで、グリッドフィルタをトップハットフィルタとみなし、さらにグリッドフィルタの特性幅について、その４乗以上に関する項を無視して、数２１で近似する。
【００９２】
【数２１】
つぎにテストフィルタを数２２で定義する。
【００９４】
【数２２】
このときテストフィルタの特性幅として、グリッドフィルタの特性幅の１乃至２倍程度の値を与えておく。
【００９６】
数２１および数２２を用いて数１３で定義されるＳＴＳ応力を整理し、グリッドフィルタの特性幅についてその４乗以上に関する項を無視すると、ＳＴＳ応力を近似するモデル式が数２３で与えられる。
【００９７】
【数２３】
数２３のモデル式は、値が未定であるモデル係数を、１つも含んでいない。
【００９９】
数１０，数１１，数１８，数１９，数２３を局所スケール間平衡仮定関係式である数１７に代入すると、数２４が得られる。
【０１００】
【数２４】
ついでＳＧＳ応力テンソルを数２５のモデル式で近似する。数２５のＳＧＳ応力テンソルを近似するモデル式は、渦粘性モデル部と、非渦粘性モデル部からなる。非渦粘性モデル部は、その値が未定であるモデル係数を含まないものとする。このとき、数２５のモデル式に含まれる、値が未定のモデル係数は渦動粘性係数のみである。
【０１０２】
【数２５】
数２５を数２４に代入して整理すると、数２６が得られる。
【０１０４】
【数２６】
数２６は数２５の唯一のモデル係数である渦動粘性係数の算出式である。数２６の右辺は既知量のみから構成されるので、数２６を用いれば、渦動粘性係数を容易に求めることができる。
【０１０６】
数２６より算出した渦動粘性係数をもってＳＧＳ応力を近似し、図３のＳ５においてＳＧＳ応力項を計算すれば、流動状況に応じた、瞬時・局所的なモデル係数の最適値をもちいて解析したことになるので、高い解析精度が期待できる。さらに本実施例によれば、グリッドスケールからサブグリッドスケールへのエネルギー伝達特性が十分に反映されるため、数値的に安定であることが期待できる。
【０１０７】
ここで本発明に係る第２の実施例について説明する。ここで数２５のＳＧＳ応力モデルの非渦粘性モデル部を数２７で与える。
【０１０８】
【数２７】
すなわちＳＧＳ応力を、渦粘性モデルのみで近似する。数２５，数２７を数２４に代入して整理すると、渦動粘性係数の算出式である数２８が得られる。
【０１１０】
【数２８】
数２８より算出した渦動粘性係数をもってＳＧＳ応力を近似し、図３のＳ５においてＳＧＳ応力項を計算すれば、第１の実施例で述べた効果に加えて、ＳＧＳ応力のモデルの、非渦粘性モデル部については計算しなくてすむため、計算時間の短縮と使用メモリ量の節約が可能となる。
【０１１２】
ここで、本発明に係る第３の実施例について説明する。ここで数２５のＳＧＳ応力モデルの非渦粘性モデル部を、数２９で与える。
【０１１３】
【数２９】
すなわち、非渦粘性モデル部を、バーディナモデルで与える。数２５，数２９を数２４に代入して整理すると、渦動粘性係数の算出式である数３０が得られる。
【０１１５】
【数３０】
数３０より算出した渦動粘性係数をもってＳＧＳ応力を近似し、図３のＳ５においてＳＧＳ応力項を計算すれば、第１の実施例で述べた効果に加えて、バーディナモデル部がＳＧＳ応力テンソル成分について高い近似精度を有することからして、さらに高い精度の解析結果が得られると期待できる。
【０１１７】
ここで、本発明に係る第４の実施例について説明する。ここで数２５のＳＧＳ応力モデルの非渦粘性モデル部を、数３１で与える。
【０１１８】
【数３１】
このとき、数３１に含まれるモデル定数の値は、あらかじめ設定しておく。数３１のモデル定数の設定については、まずＳＧＳ応力を数３２に示すモデルで与え、これに対して従来の技術で用いられてきたリリーの最小自乗法を適用し、算出する。
【０１２０】
【数３２】
したがって数３１を数２５に代入して得られるＳＧＳ応力のモデル式では、値が未定のモデル係数はやはり、渦動粘性係数のみである。渦動粘性係数の算出式は、数２５，数３１を数２４に代入して整理することにより得られ、数３３で表される。
【０１２２】
【数３３】
本実施例においても、第３の実施例と同様の効果が期待される。
【０１２４】
ここで、本発明に係る第５の実施例について説明する。ここで数２５のＳＧＳ応力モデルの非渦粘性モデル部を、数３４で与える。
【０１２５】
【数３４】
数３４は、ＳＧＳ応力をテイラー展開を用いて近似したものである。数２５，数３４を数２４に代入して整理することにより、本実施例での唯一のモデル係数である、渦動粘性係数の算出式が数３５のように得られる。
【０１２７】
【数３５】
本実施例においても、第３の実施例と同様の効果が期待される。
【０１２９】
ここで本発明に係る第６の実施例について説明する。ここで数２５のＳＧＳ応力モデルの渦動粘性係数を数４で与える。すなわち数２５のＳＧＳ応力モデルの渦粘性モデル部を、スマゴリンスキーモデルで与える。このとき、値が未定のモデル係数は、モデル定数Ｃのみである。数４，数２５を数２４に代入して整理すると、モデル定数Ｃの算出式である数３６が得られる。
【０１３０】
【数３６】
数３６より算出したモデル定数ＣをもってＳＧＳ応力を近似し、図３のＳ５においてＳＧＳ応力項を計算すれば、第１の実施例で述べたのと同様の効果が期待できる。また、本実施例において、数２５で与えられるＳＧＳ応力モデルの非渦粘性モデル部については、第２，第３，第４，第５の実施例で述べたのと同様のものが適用できる。
【０１３２】
ここで本発明に係る第７の実施例について説明する。ここではまず、グリッドフィルタをトップハットフィルタとみなし、グリッドフィルタの特性幅の４乗に関する項まで考慮して、数３７で近似する。
【０１３３】
【数３７】
つぎにテストフィルタを数２２で与え、２次中心差分法を用いて近似する。このとき、数２２の中の２階微分項は、打切り誤差について４次の項まで考慮すると、数３８のように表せる。
【０１３５】
【数３８】
すなわち、数２２を２次中心差分法で離散化すると、グリッドフィルタの特性幅の４乗に関する打ち切り誤差項があらわれるため、テストフィルタとしては、実際には数３９のフィルタを用いたことになる。
【０１３７】
【数３９】
数３７および数３９を用いて数１３で定義されるＳＴＳ応力を整理し、グリッドフィルタの特性幅についてその４乗以下に関する項までを考慮すると、ＳＴＳ応力を近似するモデル式が数４０で与えられる。
【０１３９】
【数４０】
数４０のモデル式は、値が未定であるモデル係数を、１つも含んでいない。
【０１４１】
ＳＧＳ応力のモデル式である数２５と、ＳＴＳ応力のモデル式である数４０を局所スケール間平衡仮定関係式である数１７に代入して整理すると、数４１が得られる。
【０１４２】
【数４１】
数４１は数２５の唯一のモデル係数である渦動粘性係数の算出式である。数４１の右辺は既知量のみから構成されるので、数４１を用いれば、渦動粘性係数を容易に求めることができる。
【０１４４】
数４１より算出した渦動粘性係数をもってＳＧＳ応力を近似し、図３のＳ５においてＳＧＳ応力項を計算すれば、流動状況に応じた、瞬時・局所的なモデル係数の最適値をもちいて解析したことになるので、高い解析精度が期待できる。さらに本実施例によれば、グリッドスケールからサブグリッドスケールへのエネルギー伝達特性が十分に反映されるため、数値的に安定であることが期待できる。さらに本実施例は、第１，第２，第３，第４，第５，第６の実施例がＳＴＳ応力をグリッドフィルタの特性幅の２乗に関する項まで考慮したのに対し、さらに加えて４乗に関する項につていも考慮しているため、ＳＴＳ応力に関して近似精度が高い。それにともないテストスケールからサブテストスケールへのエネルギー伝達量に関する近似精度も高くなるので、第１から第６の実施例に比べてさらに高精度の解析結果が期待できる。
【０１４５】
次に、各々のデータフローを説明する。図４は、ＳＧＳ応力テンソルの渦粘性モデル部における、渦動粘性係数を算出する際のデータフローを示すものである。ＧＳ成分の速度ベクトルデータ１１は、すでに求まっているものとする。１２は、ＧＳ成分の歪率テンソルの算出部であり、ＧＳ成分の速度ベクトル１１を用いてＧＳ成分の歪率テンソルを算出する。１３では、運動エネルギーのＧＳ成分の粘性散逸量を、１２で求めたＧＳ成分の歪率テンソルを用いて算出する。１４は、ＳＧＳ応力テンソルの非渦粘性モデル部の算出部であり、ＧＳ成分の速度ベクトル１１を用いて、例えば数２７，数２９，数３１，数３４等の近似式により算出する。１５では、ＴＳ成分の速度ベクトルをＧＳ成分の速度ベクトル１１を用いて算出する。算出に用いるテストフィルタには、例えば数２２，数３９等を適用する。１６はＴＳ成分の歪率テンソルの算出部であり、１５で求めたＴＳ成分の速度ベクトルを用いて算出する。１７では、運動エネルギーのＴＳ成分の粘性散逸量を、１６で求めたＴＳ成分の歪率テンソルを用いて算出する。１８はＳＴＳ応力テンソルの算出部で、１１のＧＳ成分の速度ベクトル、１５のＴＳ成分の速度ベクトル等を用い、例えば数２３，数４０等の近似式によってＳＴＳ応力テンソルの近似値を算出する。１９では渦動粘性係数を算出する。
【０１４６】
算出に際しては、１２のＧＳ成分の歪率テンソル、１６のＴＳ成分の歪率テンソル、１３の運動エネルギーのＧＳ成分の粘性散逸量、１７の運動エネルギーのＴＳ成分の粘性散逸量、１４のＳＧＳ応力の非渦粘性モデル部、１８のＳＴＳ応力テンソルを用いる。なお、ここで求めた渦動粘性係数は、ＳＧＳ応力テンソルの近似式である数２あるいは数２５を介し、フローチャート図３のステップＳ５において、ＳＧＳ応力項の計算で用いられる。
【０１４７】
また、図５はＳＧＳ応力テンソルの渦粘性モデル部における、渦動粘性係数を算出する際のデータフローを示すものであり、ＴＳ成分の歪率テンソル算出部１６′以外は、図４と同様である。１６′はＴＳ成分の歪率テンソルの算出部であり、１２のＧＳ成分の歪率テンソルに対し、テストフィルタ操作を行うことにより算出する（図４ではＴＳ成分の速度ベクトル１５を用いて算出）。
【０１４８】
また、図６は、ＳＧＳ応力テンソルの渦粘性モデル部における、渦動粘性係数に関するモデル定数を算出する際のデータフローを示すものである。基本的には図４のデータフローと同様である。両者の違いは、図４が渦動粘性係数を求めるデータフローであったのに対し、図６が渦粘性モデルとしてスマゴリンスキーモデルを用いる際（数４）の、モデル定数Ｃを算出するデータフローである点である。２０では渦動粘性係数に関するモデル定数を算出する。算出に際しては、１２のＧＳ成分の歪率テンソル、１６のＴＳ成分の歪率テンソル、１３の運動エネルギーのＧＳ成分の粘性散逸量、１７の運動エネルギーのＴＳ成分の粘性散逸量、１４のＳＧＳ応力の非渦粘性モデル部、１８のＳＴＳ応力テンソルを用いる。
【０１４９】
また、図７は、ＴＳ成分の歪率テンソル算出部１６′以外は、図６と同様である。
【０１５０】
以上のように、本実施例では、乱流の解析方法の一つであるＬＥＳにおいて、ＳＧＳ応力のモデル式に含まれるモデル係数の局所的な値を、瞬時の流動状況に応じて自動的に算出するため、モデル係数のチューニングを必要とせずに、高い精度の解析結果が得られるという効果を奏する。
【０１５１】
ここで、この乱流解析方法によるプログラムを作成し、このプログラム若しくはこのプログラムを記録した記録媒体でコンピュータによって容易に乱流解析を実施できる。
【０１５２】
またＬＥＳのＳＧＳ応力のモデル係数を算出するに際し、本実施例では波数空間上でのエネルギー伝達特性のモデル式である、局所スケール間平衡仮定関係式をもちいる。さらに本実施例では特に、グリッドスケールからサブグリッドスケールへのエネルギー伝達量の評価に、その定義式である数１１の関係を用いる。そのため本実施例によれば、ＳＧＳ応力のモデル係数の算出に際し、エネルギー伝達特性が十分に考慮されるから、解析が数値的に安定化するという効果が得られる。
【０１５３】
更に、テストスケールからサブテストスケールへのエネルギー伝達量を、その定義式である数１９で与えるので、高い解析精度が得られるという効果がある。
【０１５４】
更に、運動エネルギーのＴＳ成分，ＧＳ成分の粘性散逸量を、その定義式である数１０，数１８で与えるので、さらに高い解析制度が得られるという効果がある。
【０１５５】
更に、テストスケールからサブテストスケールへのエネルギー伝達量，運動エネルギーのＴＳ成分およびＧＳ成分の粘性散逸量に関し、物理的に正しい定義式を用いたことにより、高い解析精度が得られる効果がある。
【０１５６】
更に、ＳＧＳ応力を、渦粘性型のモデル式を含んだモデル式で与えるものである。渦粘性型のモデルは計算の安定化に寄与するから、本実施例によれば、数値的に極めて安定な乱流解析が可能となるという効果が得られる。
【０１５７】
【発明の効果】
本発明によると、乱流の解析方法の一つであるＬＥＳにおいて、ＳＧＳ応力のモデル式中に含まれるモデル係数もしくはモデル定数を、自動的に算出し、数値的に安定な乱流の解析方法，乱流の解析システム並びに乱流解析プログラムを提供することができるという効果を奏する。
【図面の簡単な説明】
【図１】波数空間上での運動エネルギー伝達の概念図。
【図２】局所スケール間平衡仮定の説明図。
【図３】ＬＥＳをもちいた流れの非定常計算のフローチャート。
【図４】渦動粘性係数の算出に関するデータフロー。
【図５】渦動粘性係数の算出に関するデータフロー。
【図６】渦動粘性係数に関するモデル定数の算出に関するデータフロー。
【図７】渦動粘性係数に関するモデル定数の算出に関するデータフロー。
【符号の説明】
１…運動エネルギーのＧＳ成分の粘性散逸量、２…運動エネルギーのＳＧＳ成分の粘性散逸量、３…グリッドスケールからサブグリッドスケールへのエネルギー伝達量、１′…運動エネルギーのＧＳ成分の粘性散逸量、２′…運動エネルギーのＳＧＳ成分の粘性散逸量、３′…グリッドスケールからサブグリッドスケールへのエネルギー伝達量、４…運動エネルギーのＴＳ成分の粘性散逸量、５…運動エネルギーのＳＴＳ成分でありかつＧＳ成分である成分の粘性散逸量、６…テストスケールからサブテストスケールへのエネルギー伝達量、７…サブテストスケールでありかつグリッドスケールである波数領域、１１…ＧＳ成分の速度ベクトルのデータ、１２…ＧＳ成分の歪率テンソルの算出部、１３…運動エネルギーのＧＳ成分の粘性散逸量の算出部、１４…ＳＧＳ応力テンソルの非渦粘性モデル部の算出部、１５…ＴＳ成分の速度ベクトルの算出部、１６…ＴＳ成分の歪率テンソルの算出部、１６′…ＴＳ成分の歪率テンソルの算出部、１７…運動エネルギーのＴＳ成分の粘性散逸量の算出部、１８…ＳＴＳ応力テンソルの算出部、１９…渦動粘性係数の算出部、２０…渦動粘性係数に関するモデル定数の算出部、Ｓ１…計算条件等の入力ステップ、Ｓ２…初期値設定ステップ、Ｓ３…解析時刻ゼロクリアステップ、Ｓ４…移流項の計算ステップ、Ｓ５…ＳＧＳ応力項・粘性項の計算ステップ、Ｓ６…圧力項の計算ステップ、Ｓ７…速度場の更新ステップ、Ｓ８…質量保存則の満足度の判定ステップ、Ｓ９…速度・圧力の修正計算ステップ、Ｓ１０…計算結果の出力ステップ、Ｓ１１…計算終了の判定ステップ、Ｓ１２…解析時刻の更新ステップ。
Claims (2)
- 解析しようとする流れ場を空間的に有限個の領域若しくは有限個の点に離散化し、流体の運動を非定常計算する乱流解析プログラムであって、
該離散化に伴い空間的には解像することができなくなった高波数域の速度の空間変動成分が流体の運動に与える影響であるサブグリッドスケール応力（ＳＧＳ応力）をモデル化したモデル式を計算する計算手順が前記計算機メモリもしくは記録媒体に保存され、
物理量の空間変動の波数成分のうち該離散化に伴い空間的には解像することができなくなった高波数域（サブグリッドスケール）の成分をサブグリッドスケール成分（ＳＧＳ成分），物理量の空間変動の波数成分のうちＳＧＳ成分ではない低波数域（グリッドスケール）の成分をグリッドスケール成分（ＧＳ成分）とし、
物理量の空間変動の高波数成分を除去する作用をもつ空間フィルタであるテストフィルタを、物理量の空間変動の前記ＧＳ成分に対して作用させ、前記ＧＳ成分のうちの高波数域の成分を除去した残りの低波数域(テストスケール)の成分をテストスケール成分（ＴＳ成分）とし、物理量の空間変動の波数成分のうち前記ＴＳ成分には含まれない高波数域
（サブテストスケール）の成分をサブテストスケール成分（ＳＴＳ成分）とする場合、
速度ベクトルの前記ＧＳ成分を算出し、前記計算メモリもしくは前記記録媒体に記憶する手段と、
前記算出された前記速度ベクトルのＧＳ成分のデータに基づき、運動エネルギーのＴＳ成分の粘性散逸量の近似値と、運動エネルギーのＧＳ成分の粘性散逸量の近似値と、テストスケールからサブテストスケールへのエネルギー伝達量の近似値とを算出する手段と、
運動エネルギーの該ＴＳ成分の粘性散逸量とテストスケールからサブテストスケールへのエネルギー伝達量との和が、運動エネルギーの該ＧＳ成分の粘性散逸量とグリッドスケールからサブグリッドスケールへのエネルギー伝達量との和に等しくなるという、前記計算メモリもしくは前記記録媒体に保存された局所スケール間平衡仮定関係式の計算手順を呼び出す手段と、
グリッドスケールからサブグリッドスケールへのエネルギー伝達量を、ＳＧＳ応力テンソルと、速度の空間変動のＧＳ成分に対する歪率テンソルとの内積に（−１）を乗じたのち、テストフィルタを作用させて得られる値であるとみなし、このとき該ＳＧＳ応力テンソルを、その値をまだ定めていないモデル係数を、１つだけ含む前記モデル式で与えるとき、該局所スケール間平衡仮定関係式を、該値をまだ定めていないモデル係数について整理することにより、該ＳＧＳ応力のモデル式に含まれる該モデル係数を算出するためのモデル係数算出式を導き、該モデル係数算出式を用いて該モデル係数を計算してグリッドスケールからサブグリッドスケールへのエネルギー伝達量を計算する手段とを備えたことを特徴とする乱流解析プログラム。 - 解析しようとする流れ場を空間的に有限個の領域若しくは有限個の点に離散化し、流体の運動を非定常計算する乱流解析プログラムを記録した記録媒体であって、
前記記録媒体に記録された、該離散化に伴い空間的には解像することができなくなった高波数域の速度の空間変動成分が流体の運動に与える影響であるサブグリッドスケール応力（ＳＧＳ応力）をモデル化するモデル式の計算手順が計算機メモリに保存され、
物理量の空間変動の波数成分のうち該離散化に伴い空間的には解像することができなくなった高波数域（サブグリッドスケール）の成分をサブグリッドスケール成分（ＳＧＳ成分），物理量の空間変動の波数成分のうちＳＧＳ成分ではない低波数域（グリッドスケール）の成分をグリッドスケール成分（ＧＳ成分）とし、
物理量の空間変動の高波数成分を除去する作用をもつ空間フィルタであるテストフィルタを、物理量の空間変動の前記ＧＳ成分に対して作用させ、前記ＧＳ成分のうちの高波数域の成分を除去した残りの低波数域(テストスケール)の成分をテストスケール成分（ＴＳ成分）とし、物理量の空間変動の波数成分のうち前記ＴＳ成分には含まれない高波数域
（サブテストスケール）の成分をサブテストスケール成分（ＳＴＳ成分）とする場合、
前記計算メモリから速度ベクトルの前記ＧＳ成分算出する手段と、
前記算出された速度ベクトルのＧＳ成分のデータに基づき、運動エネルギーのＴＳ成分の粘性散逸量の近似値と、運動エネルギーのＧＳ成分の粘性散逸量の近似値と、テストスケールからサブテストスケールへのエネルギー伝達量の近似値とを算出する手段と、
運動エネルギーの該ＴＳ成分の粘性散逸量とテストスケールからサブテストスケールへのエネルギー伝達量との和が、運動エネルギーの該ＧＳ成分の粘性散逸量とグリッドスケールからサブグリッドスケールへのエネルギー伝達量との和に等しくなるという、前記計算メモリに保存された局所スケール間平衡仮定関係式の計算手順を呼び出す手段と、
グリッドスケールからサブグリッドスケールへのエネルギー伝達量を、ＳＧＳ応力テンソルと、速度の空間変動のＧＳ成分に対する歪率テンソルとの内積に（−１）を乗じたのち、テストフィルタを作用させて得られる値であるとみなし、このとき該ＳＧＳ応力テンソルを、その値をまだ定めていないモデル係数を、１つだけ含む前記モデル式で与えるとき、該局所スケール間平衡仮定関係式を、該値をまだ定めていないモデル係数について整理することにより、該ＳＧＳ応力のモデル式に含まれる該モデル係数を算出するためのモデル係数算出式を導き、該モデル係数算出式を用いて該モデル係数を計算してグリッドスケールからサブグリッドスケールへのエネルギー伝達量を計算する手段とを備えたことを特徴とする乱流解析プログラムを記録した記録媒体。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2001340119A JP3882584B2 (ja) | 2001-11-06 | 2001-11-06 | 乱流の解析システム並びに乱流解析プログラム |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2001340119A JP3882584B2 (ja) | 2001-11-06 | 2001-11-06 | 乱流の解析システム並びに乱流解析プログラム |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JP2003141181A JP2003141181A (ja) | 2003-05-16 |
| JP3882584B2 true JP3882584B2 (ja) | 2007-02-21 |
Family
ID=19154368
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2001340119A Expired - Fee Related JP3882584B2 (ja) | 2001-11-06 | 2001-11-06 | 乱流の解析システム並びに乱流解析プログラム |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JP3882584B2 (ja) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107676273A (zh) * | 2017-09-25 | 2018-02-09 | 江苏大学 | 一种低振高效离心泵叶轮数值优化方法 |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN113111610B (zh) * | 2021-05-10 | 2022-10-14 | 中国空气动力研究与发展中心计算空气动力研究所 | 一种亚格子尺度模型建立方法 |
-
2001
- 2001-11-06 JP JP2001340119A patent/JP3882584B2/ja not_active Expired - Fee Related
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107676273A (zh) * | 2017-09-25 | 2018-02-09 | 江苏大学 | 一种低振高效离心泵叶轮数值优化方法 |
| CN107676273B (zh) * | 2017-09-25 | 2019-04-30 | 江苏大学 | 一种低振高效离心泵叶轮数值优化方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| JP2003141181A (ja) | 2003-05-16 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Frontiere et al. | CRKSPH–a conservative reproducing kernel smoothed particle hydrodynamics scheme | |
| Huerta et al. | A simple shock‐capturing technique for high‐order discontinuous Galerkin methods | |
| Oñate et al. | Stabilization techniques for finite element analysis of convection-diffusion problems | |
| Vishnampet et al. | A practical discrete-adjoint method for high-fidelity compressible turbulence simulations | |
| El-Kafafy et al. | Advanced frequency-domain modal analysis for dealing with measurement noise and parameter uncertainty | |
| CN109858158A (zh) | 一种计算流体力学模拟的参数配置方法及系统 | |
| Ghiasi et al. | Modal explicit filtering for large eddy simulation in discontinuous spectral element method | |
| Eguchi et al. | Numerical pressure retrieval from velocity measurement of a turbulent tornado-like vortex | |
| Van Renterghem et al. | Prediction-step staggered-in-time FDTD: an efficient numerical scheme to solve the linearised equations of fluid dynamics in outdoor sound propagation | |
| Guy et al. | A multigrid method for a model of the implicit immersed boundary equations | |
| JP3882584B2 (ja) | 乱流の解析システム並びに乱流解析プログラム | |
| Mansur et al. | Explicit time-domain approaches based on numerical Green’s functions computed by finite differences–the ExGA family | |
| Simarro et al. | A semi‐implicit non‐hydrostatic dynamical kernel using finite elements in the vertical discretization | |
| Lorentzon et al. | On stability and relaxation techniques for partitioned fluid‐structure interaction simulations | |
| KR20120132101A (ko) | 비선형 감쇠 시스템에서의 비조화적이고 주기적인 가진력 추정방법 | |
| Aricò et al. | MAST-RT0 solution of the incompressible Navier–Stokes equations in 3D complex domains | |
| Codina | A nodal‐based implementation of a stabilized finite element method for incompressible flow problems | |
| Tasri | Accuracy of cell centres to vertices interpolation for unstructured mesh finite volume solver | |
| Pahlevani et al. | Parameter Sensitivity of an eddy viscosity model: analysis, computation and its application to quantifying model reliability | |
| Moura et al. | Eigenanalysis of gradient-jump penalty (gjp) stabilisation for cg | |
| Lai et al. | A Discontinuous Galerkin Method for Two‐Dimensional Shock Wave Modeling | |
| JP7395456B2 (ja) | シミュレーション装置、及びプログラム | |
| Figueroa | Improvement of nested Cartesian finite difference grid solvers | |
| Parida et al. | Regularized numerical integration of multibody dynamics with the generalized α method | |
| Chola | SPH consistent with explicit LES |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20040119 |
|
| RD01 | Notification of change of attorney |
Free format text: JAPANESE INTERMEDIATE CODE: A7421 Effective date: 20060419 |
|
| A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20060725 |
|
| A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20060925 |
|
| 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: 20061024 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20061106 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20101124 Year of fee payment: 4 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20101124 Year of fee payment: 4 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20111124 Year of fee payment: 5 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20111124 Year of fee payment: 5 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20121124 Year of fee payment: 6 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20121124 Year of fee payment: 6 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20131124 Year of fee payment: 7 |
|
| S111 | Request for change of ownership or part of ownership |
Free format text: JAPANESE INTERMEDIATE CODE: R313111 |
|
| R350 | Written notification of registration of transfer |
Free format text: JAPANESE INTERMEDIATE CODE: R350 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| LAPS | Cancellation because of no payment of annual fees |