(57)【要約】【課題】輸送方程式における有限差分数値計算法に関し、安定性を保証する３次精度計算スキームの解析的構成法【解決手段】  輸送方程式における移流項を3次精度の有限差分近似する際に、1つのフリーパラメータεを導入する。そして安定性を保証する十分条件である差分係数の正値条件を満たし、かつ同時に安定性を損なわない範囲で、輸送速度などの局所情報量に基づいて、打切り誤差を最小とする手続きを踏んで、εの最適値を決定する。このように計算スキームを構成する際に導入したフリーパラメータεを、安定性と精度の観点からその最適値を、輸送速度、拡散係数、時間・空間メッシュ幅の局所情報量の関数として、解析的に設定することにより、上記課題を解決するもので、この点に本発明の特徴がある。【選択図】図１
2009187431.tif 000002
【特許請求の範囲】
【請求項１】
１次元、２次元、３次元移流・拡散方程式の有限差分数値計算法であって、その移流項に３次精度多項式展開差分法を適用する場合に、数値的安定性を保証する計算スキームの構成方法
【請求項２】
１次元、２次元、３次元移流・拡散方程式にソース項を加えた輸送方程式の有限差分数値計算法であって、請求項１に記載の数値的安定性を保証する計算スキームの構成方法
【請求項３】
１次元、２次元、３次元移流・拡散方程式及びそれらにソース項を加えた輸送方程式の拡散係数に対し、流体の分子粘性係数の他に、乱流効果を取り込んだ乱流粘性などを含む実効的な拡散係数を付与する場合の有限差分数値計算法であって、請求項１に記載の数値的安定性を保証する計算スキームの構成方法
【請求項４】
請求項１、請求項２、請求項３で対象とする輸送方程式は、非保存系表示、保存系表示のいずれに対しても、適用される。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、高次精度輸送数値計算の技術分野である。
【背景技術】
【０００２】
移流・拡散方程式、数１に代表される輸送方程式の数値計算の分野において、河村スキーム（非特許文献１）は、計算精度が高次の３次精度である上に、シンプルな多項式展開法に基づいて式が構成されていることから、従来から多くの理工学分野で多用されている計算技術である。
【非特許文献１】：T. Kawamura and K. Kuwabara.，”Computation of High Reynolds Number Flow around aCircular Cylinder with Surface Roughness”，AIAA-84-0340，(1984).
【０００３】
しかしながら、河村スキームを含む従来の3次精度数値計算技術では、数２に示す係数aに見られるように、有限差分方程式の係数（これを差分係数と称する）が負となり、特に輸送速度や被輸送量の空間的勾配が大きい場合に、図２の移流・拡散方程式に示す様に数値解が不安定になる。この不安定現象は元の輸送方程式の真の解には無い非物理的な数値振動である。
【０００４】
【数１】
2009187431.tif 000003
【０００５】
上式で、fは被輸送量、xは空間座標、uは輸送速度（以下の説明では正を仮定）、νは拡散係数を表す。上式の左辺第２項を移流項、右辺を拡散項と称する。
【０００６】
移流・拡散方程式、数１の時間項をオイラー陽解法で、移流項を河村スキームで、拡散項を４次精度中心差分スキームで離散近似した時の、被輸送量fの新時刻ステップ（n+1）における計算式は、数２の河村スキームの有限差分方程式と差分係数で与えられる。
【０００７】
【数２】
2009187431.tif 000004
【０００８】
ここで、C、Dは次式で定義される無次元数である。
【０００９】
【数３】
2009187431.tif 000005
【００１０】
数３で、上添え字nは時刻ステップ番号、下添え字iは計算対象の空間を計算格子に分割した時のメッシュ番号である。
【００１１】
図２は移流・拡散方程式、数１を従来技術の河村スキームで計算した時の数値解である。
計算条件：空間メッシュ総数=100、C=0.17、D=0.25、時刻ステップ数300
【発明の開示】
【発明が解決しようとする課題】
【００１２】
従来の３次精度多項式展開差分技術で輸送方程式を計算する場合に、解消されない数値的不安定性問題。
【００１３】
本発明は、従来の高次精度輸送計算法における不安定性の問題を解決しようとするもので、精度は従来同様３次精度を維持しつつ、同時に安定性を確実に保証するものである。
【００１４】
従来、安定性が損なわれる場所では、計算スキームの精度を落とし、そのことによって必然的に付加される数値誤差（数値拡散という）によって解の不安定性を回避する手法が一般的に採用されていた。
【００１５】
しかしこの回避法では、衝撃波のような急峻な勾配を有する問題の場合に勾配の鈍化現象が発生し、数値解の質が低下してしまう。
【００１６】
そこで本発明は、いかなる場所においても計算スキームの精度を落とすことなく３次精度を維持し、安定性を損なわない範囲で打切り誤差を最小とするような最適なパラメータを解析的に決定するものである。
【課題を解決するための手段】
【００１７】
輸送方程式における移流項に関し、３次精度の計算スキームを表現する際に、１つのフリーパラメータεを導入する。安定性を保証する十分条件である差分係数の正値条件を満たし、かつ輸送速度、拡散係数、時間・空間メッシュ幅の局所情報量に基づいて、安定性を損なわない範囲で打切り誤差を最小とする手続きを踏んで、εを決定する。
【００１８】
このように計算スキームを表現する際に導入したフリーパラメータεを、安定性と精度の観点から最適値に設定することにより、上記課題を解決するものである。
【００１９】
移流項に対し、フリーパラメータεを導入した３次精度計算スキームの一般式を数３に示す。
【００２０】
【数３】
2009187431.tif 000006
【００２１】
上式は、ε＝0.5の時河村スキーム、ε=0の時UTOPIAスキームに一致する。
上式の差分係数を数４に示す。
【００２２】
【数４】
2009187431.tif 000007
【発明の効果】
【００２３】
数６、数７で決まるパラメータεの最適値を用いて、図２と同一条件で計算した結果を、図１に示す。最適化スキームによって、安定な数値解を得ることができる。計算時間は僅かに増加する傾向であるが、実計算上、問題になることはない。
【発明を実施するための最良の形態】
【００２４】
数１の輸送方程式の時間項をオイラー法で、移流項を、数３、数４の式で、拡散項を４次精度中心差分式で近似して得られる差分方程式と、その差分係数a、b、c、d、eを数５に示す
【００２５】
数値解の安定性を保証する十分条件は、差分係数a、b、c、d、eの正値条件から得られる。この正値条件から得られるクーラン数Cと拡散数Dに対する許容条件を数６に、更に打切り誤差TEの係数の絶対値が最小となる条件より決定されるεの最適値を数７に示す。
【００２６】
数７に示すように、フリーパラメータεの最適値を、時間メッシュ幅Δt、空間メッシュ幅Δx、輸送速度u、拡散係数νの局所情報量の関数として解析的に決定したことに、本発明の特徴がある。このスキームを、最適化スキームと称する。
【００２７】
多次元輸送方程式の場合も同様にして、最適スキームを決定することができる。２次元輸送方程式の場合のC、Dの許容範囲を数８とパラメータεの最適値を数９に、３次元輸送方程式のC、Dの許容範囲を数１０とパラメータεの最適値を数１１に示す。
【００２８】
【数５】
2009187431.tif 000008
【００２９】
【数６】
2009187431.tif 000009
【００３０】
【数７】
2009187431.tif 000010
【００３１】
【数８】
2009187431.tif 000011
【００３２】
【数９】
2009187431.tif 000012
【００３３】
【数１０】
2009187431.tif 000013
【００３４】
【数１１】
2009187431.tif 000014
【図面の簡単な説明】
【００３５】
【図１】移流・拡散方程式、数１を本発明の最適化スキームで計算したときの数値解を示した図である。
【図２】移流・拡散方程式、数１を従来技術の河村スキームで計算したときの数値解である。計算条件：空間メッシュ総数＝１００、C＝０．１７、D＝０．２５、時刻ステップ数３００
【図１】
2009187431.tif 000015
【図２】
2009187431.tif 000016