(57)【要約】      （修正有）【課題】シミュレーションに必要なメモリ及び計算時間を大幅に低減できる物理量シミュレーション方法を提供する。【解決手段】解析対象領域を等サイズの矩形ブロックに分割し、互いに第１間隔で離間する基本グリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第１格子点のそれぞれに対応した物理量をＣＩＰ法により演算し、前記基本グリッドと一部のグリッドが重複すると共に、互いに前記第１間隔より短い第２間隔で離間するサブグリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第２格子点のそれぞれに対応した物理量をＣＩＰ法により演算し、前記矩形ブロックに含まれる格子点に対応する物理量が所定の閾値以上でないと判断された場合には、前記第１演算ステップによる演算を実行する共に、物理量が所定の閾値以上であると判断された場合には、前記第２演算ステップによる演算を実行する。【選択図】図６
2013214299.tif 000002
【特許請求の範囲】
【請求項１】
解析対象領域における物理量をシミュレーションする物理量シミュレーション方法であって、
解析対象領域を等サイズの矩形ブロックに分割する分割ステップと、
前記解析対象領域において、互いに第１間隔で離間する基本グリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第１格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第１演算ステップと、
前記解析対象領域において、前記基本グリッドと一部のグリッドが重複すると共に、互いに前記第１間隔より短い第２間隔で離間するサブグリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第２格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第２演算ステップと、
前記矩形ブロックに含まれる格子点に対応する物理量が所定の閾値以上であるか否かを判定する判定ステップと、
前記判定ステップにより物理量が所定の閾値以上でないと判断された場合には、前記第１演算ステップによる演算を実行する共に、前記判定ステップにより物理量が所定の閾値以上であると判断された場合には、前記第２演算ステップによる演算を実行する演算切り換えステップと、を有することを特徴とする物理量シミュレーション方法。
【請求項２】
解析対象領域における物理量をシミュレーションする物理量シミュレーション方法であって、
解析対象領域を等サイズの矩形ブロックに分割する分割ステップと、
前記解析対象領域において、互いに第１間隔で離間する基本グリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第１格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第１演算ステップと、
前記解析対象領域において、前記基本グリッドと一部のグリッドが重複すると共に、互いに前記第１間隔より短い第２間隔で離間するサブグリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第２格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第２演算ステップと、
前記矩形ブロックに含まれる格子点に対応する物理量が所定の閾値以上であるか否かを判定する判定ステップと、
前記判定ステップにより物理量が所定の閾値以上でないと判断された場合には、前記矩形ブロックに含まれる格子点について前記第１演算ステップによる演算を実行する共に、前記判定ステップにより物理量が所定の閾値以上であると判断された場合には、前記矩形ブロックに含まれる格子点と、前記矩形ブロックと隣接する矩形ブロックに含まれる格子点について前記第２演算ステップによる演算を実行する演算切り換えステップと、を有し、前記演算切り換えステップは、前記矩形ブロックの１辺の長さと、物理量が解析対象領域を伝搬する速度とに応じて決まる時間毎に実行されることを特徴とする物理量シミュレーション方法。
【請求項３】
解析対象領域における特定方向に伝搬する波動に係る物理量をシミュレーションする物理量シミュレーション方法であって、
解析対象領域を等サイズの矩形ブロックに分割する分割ステップと、
前記解析対象領域において、互いに第１間隔で離間する基本グリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第１格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第１演算ステップと、
前記解析対象領域において、前記基本グリッドと一部のグリッドが重複すると共に、互いに前記第１間隔より短い第２間隔で離間するサブグリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第２格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第２演算ステップと、
前記矩形ブロックに含まれる格子点に対応する、ＣＩＰ法により物理量を演算する際に利用される、特定方向に伝搬する波動の成分を表す特性値が所定の閾値以上であるか否かを判定する判定ステップと、
前記判定ステップにより特性値が所定の閾値以上でないと判断された場合には、前記矩形ブロックに含まれる格子点について前記第１演算ステップによる演算を実行する共に、前記判定ステップにより特性値が所定の閾値以上であると判断された場合には、前記矩形ブロックに含まれる格子点と、前記特定方向で前記矩形ブロックと隣接する矩形ブロックに含まれる格子点について前記第２演算ステップによる演算を実行する演算切り換えステップと、を有し、
前記演算切り換えステップは、前記矩形ブロックの１辺の長さと、物理量が解析対象領域を伝搬する速度とに応じて決まる時間毎に実行されることを特徴とする物理量シミュレーション方法。
【請求項４】
前記物理量が音圧であることを特徴とする請求項１乃至請求項３のいずれか１項に記載の物理量シミュレーション方法。
【請求項５】
前記物理量が粒子速度であることを特徴とする請求項１乃至請求項３のいずれか１項に記載の物理量シミュレーション方法。
【請求項６】
請求項１乃至請求項５のいずれかに記載の方法を実行することを特徴とする物理量シミュレーションシステム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、コンピューター上で音波や電磁波の時間に伴う伝搬をシミュレーションするために好適な物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムに関する。
【背景技術】
【０００２】
従来より、コンピューター上で音波や電磁波の時間に伴う伝搬を波動理論に基づきシミュレーションする手法として、ＣＩＰ（Ｃｏｎｓｔｒａｉｎｅｄ  Ｉｎｔｅｒｐｏｌａｔｉｏｎ  Ｐｒｏｆｉｌｅ）法が提案されている（例えば、非特許文献１）。
【０００３】
上記のようなＣＩＰ法では、波動の伝搬計算に場の物理量（例えば、音波伝搬では音圧と粒子速度）に加えてその方向微分値も用いる。このため、グリッドで定義される格子点に対応する場の物理量とその方向微分値をメモリに保持する必要があり、ＣＩＰ法によるシミュレーションではコンピューターの使用メモリが多い。
【０００４】
また、ＣＩＰ法によるシミュレーションでは、波動の伝搬過程で実際の物理現象以上に振幅が減衰する「数値散逸」の問題があり、伝搬距離の大きい大規模空間や屋外の波動伝搬解析へ適用する場合、その影響は無視できない。グリッド数を増やして空間を細かく離散化すれば数値散逸の影響は低減できるが、グリッド数の増加は、上述のように使用メモリの増大を招き、計算時間も増大する。
【０００５】
そこで、数値散逸の課題を解決する手段の一つとして、同じく時間領域解析法であるＦＤＴＤ（Ｆｉｎｉｔｅ  Ｄｉｆｆｅｒｎｃｅ  Ｔｉｍｅ  Ｄｏｍａｉｎ）法においても用いられる、サブグリッド手法を適用したＣＩＰ法（以下、サブグリッドＣＩＰ法）を用いたシミュレーションが提案されている（例えば、非特許文献２）。このサブグリッドＣＩＰ法では、解析対象領域のうち高い精度が必要な領域を細かく離散化（サブグリッド）し、その他の領域は粗く離散化（基本グリッド）する方法で全体のグリッド数を抑制する。
【非特許文献１】「Ｃ型ＣＩＰ法を用いた音場解析に関する検討」、信学技報、Ｖｏｌ．１０６，Ｎｏ．４８１，ｐｐ．１７−２２，ｓ２００６−９８，２００６．
【非特許文献２】「ＣＩＰ法による音波伝搬シミュレーションにおけるサブグリッド・テクニックの精度比較」，日本音響学会２０１１年秋季研究発表会講演論文集，ｐｐ．１４９７−１５００，２０１１．
【発明の概要】
【発明が解決しようとする課題】
【０００６】
しかしながら、従来のサブグリッドＣＩＰ法を用いたシミュレーションでは、伝搬する音波の波面が時間の経過と共に移動して、解析対象領域が細かく離散化されているサブグリッド領域の外に出た場合は、計算精度が低下する、という問題があった。そこで、解析対象領域が細かく離散化されているサブグリッド領域を、解析対象領域内において広く取るようにすれば、広い範囲で計算精度を保つこと可能となる。しかしながら、このようにすると、全体のグリッド数の増加に伴い格子点が増えるので、シミュレーションに必要なメモリ及び計算時間が増大する、という新たな問題が生じることとなる。
【課題を解決するための手段】
【０００７】
この発明は、上記課題を解決するものであって、請求項１に係る発明は、解析対象領域における物理量をシミュレーションする物理量シミュレーション方法であって、解析対象
領域を等サイズの矩形ブロックに分割する分割ステップと、前記解析対象領域において、互いに第１間隔で離間する基本グリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第１格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第１演算ステップと、前記解析対象領域において、前記基本グリッドと一部のグリッドが重複すると共に、互いに前記第１間隔より短い第２間隔で離間するサブグリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第２格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第２演算ステップと、前記矩形ブロックに含まれる格子点に対応する物理量が所定の閾値以上であるか否かを判定する判定ステップと、前記判定ステップにより物理量が所定の閾値以上でないと判断された場合には、前記第１演算ステップによる演算を実行する共に、前記判定ステップにより物理量が所定の閾値以上であると判断された場合には、前記第２演算ステップによる演算を実行する演算切り換えステップと、を有することを特徴とする
また、請求項２に係る発明は、解析対象領域における物理量をシミュレーションする物理量シミュレーション方法であって、解析対象領域を等サイズの矩形ブロックに分割する分割ステップと、前記解析対象領域において、互いに第１間隔で離間する基本グリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第１格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第１演算ステップと、前記解析対象領域において、前記基本グリッドと一部のグリッドが重複すると共に、互いに前記第１間隔より短い第２間隔で離間するサブグリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第２格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第２演算ステップと、前記矩形ブロックに含まれる格子点に対応する物理量が所定の閾値以上であるか否かを判定する判定ステップと、前記判定ステップにより物理量が所定の閾値以上でないと判断された場合には、前記矩形ブロックに含まれる格子点について前記第１演算ステップによる演算を実行する共に、前記判定ステップにより物理量が所定の閾値以上であると判断された場合には、前記矩形ブロックに含まれる格子点と、前記矩形ブロックと隣接する矩形ブロックに含まれる格子点について前記第２演算ステップによる演算を実行する演算切り換えステップと、を有し、前記演算切り換えステップは、前記矩形ブロックの１辺の長さと、物理量が解析対象領域を伝搬する速度とに応じて決まる時間毎に実行されることを特徴とする。
【０００８】
また、請求項３に係る発明は、解析対象領域における特定方向に伝搬する波動に係る物理量をシミュレーションする物理量シミュレーション方法であって、解析対象領域を等サイズの矩形ブロックに分割する分割ステップと、前記解析対象領域において、互いに第１間隔で離間する基本グリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第１格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第１演算ステップと、前記解析対象領域において、前記基本グリッドと一部のグリッドが重複すると共に、互いに前記第１間隔より短い第２間隔で離間するサブグリッドに基づいて定義され、各々の前記矩形ブロックに等しい数含まれる第２格子点のそれぞれに対応した物理量をＣＩＰ法により演算する第２演算ステップと、前記矩形ブロックに含まれる格子点に対応する、ＣＩＰ法により物理量を演算する際に利用される、特定方向に伝搬する波動の成分を表す特性値が所定の閾値以上であるか否かを判定する判定ステップと、前記判定ステップにより特性値が所定の閾値以上でないと判断された場合には、前記矩形ブロックに含まれる格子点について前記第１演算ステップによる演算を実行する共に、前記判定ステップにより特性値が所定の閾値以上であると判断された場合には、前記矩形ブロックに含まれる格子点と、前記特定方向で前記矩形ブロックと隣接する矩形ブロックに含まれる格子点について前記第２演算ステップによる演算を実行する演算切り換えステップと、を有し、前記演算切り換えステップは、前記矩形ブロックの１辺の長さと、物理量が解析対象領域を伝搬する速度とに応じて決まる時間毎に実行されることを特徴とする。
【０００９】
また、請求項４に係る発明は、請求項１乃至請求項３のいずれか１項に記載の物理量シ
ミュレーション方法において、前記物理量が音圧であることを特徴とする。
【００１０】
また、請求項５に係る発明は、請求項１乃至請求項３のいずれか１項に記載の物理量シミュレーション方法において、前記物理量が粒子速度であることを特徴とする。
【００１１】
また、請求項６に係る発明は、請求項１乃至請求項５のいずれかに記載の方法を実行することを特徴とする物理量シミュレーションシステムである。
【発明の効果】
【００１２】
本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、物理量が所定の閾値以上となる主要な波面の付近のみサブグリッドに基づく第２演算を実行することで、全体の格子点数を抑えてシミュレーションに必要なメモリ及び計算時間を大幅に低減できる。
【００１３】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、波面の移動に合わせてサブグリッドに基づく演算を行う領域を設定するので、領域全体をサブグリッドと同じ離散化幅で細かく離散化した場合と同等程度の計算精度が得られる。
【００１４】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、解析対象領域を矩形のブロックに分割することで、コンピューター内に確保したデータ配列とブロックの形が一致し、ブロック内及び境界部分の伝搬計算のアルゴリズムが簡単になる。
【００１５】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、解析対象領域を等サイズの矩形ブロックに分割することで、ブロック境界において隣接するブロックが１体１で対応するので、サブグリッドに基づく演算への切り換え時及び境界部分の伝搬計算に必要な隣接ブロックのデータを参照するアルゴリズムが簡単になる。
【００１６】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、ブロックサイズを適切に設定することにより、様々な計算条件に柔軟に対応できる。例えば、計算時間はかかっても使用メモリを節約したい場合は、小さいブロックサイズを設定する。サブグリッド領域を最小限に設定して必要メモリを少なくすることができる。この場合、ブロック数が多いためグリッド変換やブロック境界の伝搬計算にかかる計算時間は増える。逆に、使用メモリより計算時間を節約したい場合は、大きいブロックサイズを設定する。波面付近のやや広い範囲にサブグリッド領域を設定するため必要なメモリは増えるが、ブロック数が少ないためグリッド変換やブロック境界の伝搬計算にかかる計算時間は減る。
【００１７】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、基本グリッドにより演算を行うか、サブグリッドにより演算を行うか、の判定と、その切り換えを行う時間ステップ間隔を大きく取れば、計算全体で判定と変換の回数を減らせ、グリッド変換処理による計算時間の増加を少なくすることができる。
【００１８】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、可変グリッドＣＩＰ法による波動伝搬シミュレーションにおいて、基本グリッドによる演算を行うかサブグリッドにより演算を行うか、の判定と、その切り換えを行う時間ステップ間隔を大きくとることができる。また、前記判定とグリッド変
換処理の回数を減らすことにより、全体の計算時間の増加を少なくすることができる。
【００１９】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムは、音波の伝搬解析以外にも、電磁波など他の波動一般の伝搬解析にも適用可能である。
【図面の簡単な説明】
【００２０】
【図１】本発明の実施形態に係る物理量シミュレーションシステムを構成するコンピューターの一例を示す図である。
【図２】解析対象領域における矩形ブロックへの分割を説明する図である。
【図３】基本グリッド、サブグリッド間の変換を説明する図である。
【図４】エルミート補間を説明するための図である。
【図５】第１演算を行うブロックと、第２演算を行うブロックとの境界部における演算方法を説明する図である。
【図６】本発明に係る物理量シミュレーション方法をコンピューターに実行させるためのフローチャート例を示す図である。
【図７】本発明に係る方法と比較例に係る方法のシミュレーション結果を示す図である。
【図８】本発明に係る方法と比較例に係る方法のシミュレーション結果を示す図である。
【図９】本発明に係る方法と比較例に係る方法のシミュレーション結果を示す図である。
【図１０】本発明に係る物理量シミュレーション方法の概念を示す図である。
【図１１】第２実施形態に係る物理量シミュレーション方法の概念を示す図である。
【図１２】第２実施形態に係る物理量シミュレーション方法をコンピューターに実行させるためのフローチャート例を示す図である。
【図１３】第３実施形態に係る物理量シミュレーション方法の概念を示す図である。
【図１４】第３実施形態に係る物理量シミュレーション方法をコンピューターに実行させるためのフローチャート例を示す図である。
【発明を実施するための形態】
【００２１】
以下、本発明の実施の形態を図面を参照しつつ説明する。図１は本発明の実施形態に係る物理量シミュレーションシステムを構成するコンピューターの一例を示す図である。図１において、１０はシステムバス、１１はＣＰＵ（Ｃｅｎｔｒａｌ  Ｐｒｏｃｅｓｓｉｎｇ  Ｕｎｉｔ）、１２はＲＡＭ（Ｒａｎｄｏｍ  Ａｃｃｅｓｓ  Ｍｅｍｏｒｙ）、１３はＲＯＭ（Ｒｅａｄ  Ｏｎｌｙ  Ｍｅｍｏｒｙ）、１４は外部情報機器との通信を司る通信制御部、１５はキーボードコントローラなどの入力制御部、１６は出力制御部、１７は外部記憶装置制御部、１８はキーボード、ポインティングデバイス、マウスなどの入力機器からなる入力部、１９は印刷装置などの出力部、２０はＨＤＤ（Ｈａｒｄ  Ｄｉｓｋ  Ｄｒｉｖｅ）等の外部記憶装置、２１はグラフィック制御部、２２はディスプレイ装置をそれぞれ示している。
【００２２】
図１において、ＣＰＵ１１は、ＲＯＭ１３内のプログラム用ＲＯＭ、或いは、大容量の外部記憶装置２０に記憶されたプログラム等に応じて、外部機器と通信することでデータを検索・取得したり、また、図形、イメージ、文字、表等が混在した出力データの処理を実行したり、更に、外部記憶装置２０に格納されているデータベースの管理を実行したり、などといった演算処理を行うものである。
【００２３】
また、ＣＰＵ１１は、システムバス１０に接続される各デバイスを統括的に制御する。ＲＯＭ１３内のプログラム用ＲＯＭあるいは外部記憶装置２０には、ＣＰＵ１１の制御用
の基本プログラムであるオペレーティングシステムプログラム（以下ＯＳ）等が記憶されている。また、ＲＯＭ１３あるいは外部記憶装置２０には出力データ処理等を行う際に使用される各種データが記憶されている。メインメモリーであるＲＡＭ１２は、ＣＰＵ１１の主メモリ、ワークエリア等として機能する。
【００２４】
入力制御部１５は、キーボードや不図示のポインティングデバイスからの入力部１８を制御する。また、出力制御部１６は、プリンタなどの出力部１９の出力制御を行う。
【００２５】
外部記憶装置制御部１７は、ブートプログラム、各種のアプリケーション、フォントデータ、ユーザーファイル、編集ファイル、プリンタドライバ等を記憶するＨＤＤ（Ｈａｒｄ  Ｄｉｓｋ  Ｄｒｉｖｅ）や、或いはフロッピーディスク（ＦＤ）等の外部記憶装置２０へのアクセスを制御する。本発明の物理量シミュレーション方法を実現するシステムプログラムは、上記のような外部記憶装置２０に記憶されている。また、グラフィック制御部２１は、ディスプレイ装置２２に表示する情報を描画処理するための構成である。
【００２６】
また、通信制御部１４は、ネットワークを介して、外部機器と通信を制御するものであり、これによりシステムが必要とするデータを、インターネットやイントラネット上の外部機器が保有するデータベースから取得したり、外部機器に情報を送信したりすることができるように構成される。
【００２７】
外部記憶装置２０には、ＣＰＵ１１の制御プログラムであるオペレーティングシステムプログラム（以下ＯＳ）以外に、本発明の物理量シミュレーションシステムをＣＰＵ１１上で動作させるシステムプログラム、及びこのシステムプログラムで用いるデータなどがインストールされ保存・記憶されている。
【００２８】
本発明の物理量シミュレーション方法を実現するシステムプログラムで利用されるデータとしては、基本的には外部記憶装置２０に保存されていることが想定されているが、場合によっては、これらのデータを、通信制御部１４を介してインターネットやイントラネット上の外部機器から取得するように構成することも可能である。また、本発明の物理量シミュレーション方法を実現するシステムプログラムで利用されるデータを、ＵＳＢメモリやＣＤ、ＤＶＤなどの各種メディアから取得するように構成することもできる。
【００２９】
次に、上記のようなシステム構成のコンピューターにより実行可能な本発明に係る物理量シミュレーション方法について、以下説明する。
【００３０】
以下の例では、上記のようなコンピューター上で、解析対象領域における音波の時間に伴う伝搬をシミュレーションする例を説明するが、本発明は電磁波など他の波動一般の伝搬解析にも適用可能である。
【００３１】
解析対象領域における音波の時間に伴う伝搬をシミュレーションする基礎となる方程式について説明する。
【００３２】
損失を無視した場合、音場の２つの支配方程式（ｅｑｕａｔｉｏｎ  ｏｆ  ｃｏｎｔｉｎｕｉｔｙ、ｅｑｕａｔｉｏｎ  ｏｆ  ｍｏｔｉｏｎ）は、以下の式（１）及び式（２）で示される。
【００３３】
【数１】
2013214299.tif 000003
【００３４】
【数２】
2013214299.tif 000004
ここで、ρは媒質の密度、Ｋは体積弾性率、ｐは音圧、ｖは粒子速度である。
【００３５】
簡単のため１次元領域での伝搬解析とすると、式（１）と式（２）から以下の式（３）及び式（４）が得られる。
【００３６】
【数３】
2013214299.tif 000005
【００３７】
【数４】
2013214299.tif 000006
また、式（３）と式（４）の和と差を計算して、
【００３８】
【数５】
2013214299.tif 000007
が得られる。ただし、Ｚは特性インピーダンス、ｃは媒質中の音速である（Ｚ＝（ρＫ）1/2、ｃ＝（Ｋ／ρ）1/2）。
【００３９】
さらにｆ±＝ｐ±Ｚｖxとおくと
【００４０】
【数６】
2013214299.tif 000008
を得る。式（６）は±ｘ方向へ速度ｃで伝搬する波ｆ±に対する移流方程式になっている。ＣＩＰ法は場の物理量に対する微分値も用いる手法であるので、微分値についても同様に定式化を行う。式（６）をｘで偏微分することにより以下の式が与えられる。
【００４１】
【数７】
2013214299.tif 000009
ただし、
【００４２】
【数８】
2013214299.tif 000010
である。これらの式により、±ｘ方向への音波伝搬を解くことができる。また、
【００４３】
【数９】
2013214299.tif 000011
で考えると、ｙ方向についてもｘ方向と同様に計算することができる。また、３次元の音波の伝搬ついては、ｚ方向に同様の計算を追加することで計算可能となる。
【００４４】
以上のような、移流方程式に基づく音波の伝搬ついては、非特許文献１及び非特許文献２に記載のＣＩＰ法によりシミュレーションを行うことが可能である。ここでは、非特許文献１及び非特許文献２に記載のＣＩＰ法に関する記載を参照して援用する。
【００４５】
次に、ＣＩＰ法で解析対象領域において音波の伝搬のシミュレーションを行う場合の前提となる事項について説明する。
【００４６】
図２は解析対象領域における矩形ブロックへの分割を説明する図である。本発明に係る物理量シミュレーション方法においては、音波の伝搬をシミュレーションする解析対象領域については、等サイズの矩形ブロックに分割する。このように、解析対象領域を矩形のブロックに分割することで、コンピューター内に確保したデータ配列とブロックの形が一致し、ブロック内及び境界部分の伝搬計算のアルゴリズムが簡単になるというメリットがある。
【００４７】
また、解析対象領域には、点線で示すようなグリッドが定義されており、本発明に係る物理量シミュレーション方法においては、このグリッドが交差する点として定義される各格子点における物理量（２次元の音波の伝搬の場合は、音圧ｐと粒子速度ｖx、ｖy）がシミュレーションされる。上記のようなグリッドは、互いに所定の第１の間隔で離間するようにされている。また、図２では、各矩形ブロックに２×２の格子点が定義されている例を示しているが、矩形ブロック中に定義される格子点の数が特にこのように限定されるわけではない。このような格子点の数は、各々の矩形ブロックで等しくなるように設定される。なお、上記のような基本となるグリッドが交差する点である格子点を第１格子点と称する。
【００４８】
本発明においては、上記のような基本グリッドに加え、互いの間隔が前記の第１間隔より短い第２の間隔で離間するサブグリッドが定義される。このサブグリッドが交差する点として定義される各格子点（第２格子点と称する）により、ブロック内を細かく離散化することが可能となり、これにより、基本グリッドに基づく場合より高精度の計算が可能と
なる。
【００４９】
本発明に係る物理量シミュレーション方法においては、ある時間のステップにおいて、ブロック内の基本グリッドに基づく格子点の物理量（例えば、音圧等。或いはその方向微分値）が設定した閾値以上となった場合、当該ブロック内を細かく離散化したサブグリッドに基づく演算に切り換え、そのブロック内では小さい離散化幅で高精度の計算を行う。サブグリッドに基づく第２格子点のそれぞれに対応した物理量をＣＩＰ法により演算する演算を第２演算と称することとする。
【００５０】
また、逆に、ある時間のステップにおいて、ブロック内のサブグリッドに基づく格子点の物理量（例えば、音圧等。或いはその方向微分値）が設定した閾値を超えない場合には、当該ブロック内を粗く離散化した基本グリッドに基づく演算に切り換え、そのブロック内では大きい離散化幅で計算を行い、メモリの使用量や計算時間を低減する。基本グリッドに基づく第１格子点のそれぞれに対応した物理量をＣＩＰ法により演算する演算を第１演算と称することとする。
【００５１】
ここで、本発明に係る物理量シミュレーション方法における基本グリッドに基づく演算、サブグリッドに基づく演算の間の切り換えについて図３を参照して説明する。
【００５２】
図３（Ａ）は基本グリッドとそれに基づく第１格子点を示している。このような基本グリッドに基づく第１格子点は、例えば、ｘ方向においてΔｘの間隔で離間し、ｙ方向においてΔｙの間隔で離間しているものとする。基本グリッドに基づく演算においては、第１格子点に対する物理量がＣＩＰ法により演算されることとなる。
【００５３】
本発明に係る物理量シミュレーション方法においては、第１格子点の物理量が前記したように所定の閾値以上となった場合には、細かく離間化したサブグリッドに基づく演算に切り換えるが、このとき、第１格子点群の物理量から、サブグリッドに基づく第２格子点群の物理量を生成する必要がある。ここでは、サブグリッドに基づく第２格子点は、例えば、ｘ方向においてΔｘsの間隔で離間し、ｙ方向においてΔｙsの間隔で離間しているものとする。
【００５４】
図３（Ｂ）は上記のような基本グリッドに基づく第１格子点における物理量から、サブグリッドに基づく第２格子点の物理量を生成する様子を示している。ここで、第１格子点における物理量から、第２格子点の物理量を生成する際には、例えば線形補間などを代表とする各種の補間法によってこれを行うことができる。本実施形態においては、３次のエルミート補間を用いた。前述したように、ＣＩＰ法においては、物理量とその方向微分値が用いられるが、３次のエルミート補間においても、物理量とその方向微分値を利用するために、本発明に好適な方法である、ということができる。
【００５５】
図４はエルミート補間を説明するための図である。図４において、位置ｘ―における物理量ｆ―と、位置ｘ+における物理量ｆ+とから、位置ｘiにおける物理量ｆiをエルミート補間により生成することを示している。
【００５６】
３次のエルミート補間によれば、位置ｘiにおける物理量ｆiと、物理量ｆiの微分値ｄ
ｘｆiは、Ｃ1乃至Ｃ8を係数として、下式（１０）、（１１）のように現すことができる
。
【００５７】
【数１０】
2013214299.tif 000012
【００５８】
【数１１】
2013214299.tif 000013
ここで、物理量として音圧ｐを扱う場合は、ｆ→ｐとして、物理量として粒子速度ｖx
を扱う場合は、ｆ→ｖxとして、物理量として粒子速度ｖyを扱う場合は、ｆ→ｖyとする
。
【００５９】
図３（Ｂ）において、第２格子点（１）における物理量は、ｘ方向両側の第１格子点の物理量からエルミート補間により値を生成する。また、第２格子点（２）における物理量は、ｘ方向両側の第１格子点の物理量からエルミート補間により値を生成する。また、第２格子点（３）における物理量は、ｙ方向両側の第１格子点の物理量からエルミート補間により値を生成する。また、第２格子点（４）における物理量は、ｙ方向両側の第１格子点の物理量からエルミート補間により値を生成する。また、第２格子点（５）における物理量は、第２格子点（１）及び第２格子点（２）における物理量から、エルミート補間により値を生成する。
【００６０】
以上、図３（Ｂ）に示すような手順により、本発明に係る物理量シミュレーション方法においては、第１格子点における物理量から、第２格子点の物理量をエルミート補間により生成する。
【００６１】
また、本発明に係る物理量シミュレーション方法においては、第２格子点の物理量が前記したように所定の閾値を下回った場合には、細かく離散化したサブグリッドに基づく第２演算から、粗く離散化した基本グリッドに基づく第１演算に切り換える。このような場合は、図３（Ｃ）に示すように、必要となる第２格子点の物理量を抽出し、その他の第２格子点の物理量は破棄してメモリを解放する。
【００６２】
ところで、隣接するブロックのグリッドのタイプが同じ場合、境界部分の離散化幅は隣接するブロックで同じなので、通常のＣＩＰ法による伝搬シミュレーションを行うことができる。
【００６３】
しかしながら、隣接ブロックのグリッドタイプが異なり、一方のブロックでは第１演算が行われ、他方のブロックでは第２演算が行われる場合、境界部分の離散化幅が異なる。このため基本グリッドに基づいて格子点の演算を行っているブロックにおいても、補間によって必要となる格子点における物理量を求めて、伝搬シミュレーションを行うようにする。
【００６４】
図５は第１演算を行うブロックと、第２演算を行うブロックとの境界部における演算方法を説明する図である。図５に示すように、第２演算を行うブロックで必要となる、白丸で示す第２格子点の物理量については、第１演算を行うブロックにおいても、エルミート補間により値を生成して、第２演算を行うブロックで利用するようにしている。
【００６５】
なお、ＣＩＰ法による２次元或いは３次元空間のシミュレーションでは、場の物理量の
方向微分値は、伝搬軸方向に加えそれ以外の軸方向への微分値も計算に加える必要がある。例えば、ｘｙｚ座標系では、ｘ方向の伝搬計算に、ｘ方向微分値に加えてｙ、ｚ方向微分値についても必要とする。また、特に、Ｃ型ＣＩＰ法では、更にｘｙ方向への２階微分値、ｘｚ方向への２階微分値も計算に必要になる。新たに生成されたグリッドにおけるこれらの伝搬軸方向以外の方向微分値については、前記したようにエルミート補間を使って求めてもよいが、計算量の少ない線形補間により求めてもよい。補間の精度はやや低下するが、一般に伝搬軸方向以外の軸方向の微分値は大きく変化しないので、全体の計算精度への影響は小さい。
【００６６】
次に、本発明に係る物理量シミュレーション方法をコンピューターに実行させる際のアルゴリズムについて説明する。図６は本発明に係る物理量シミュレーション方法をコンピューターに実行させるためのフローチャート例を示す図である。
【００６７】
ここで、フローチャート中で用いられる変数・定数について説明する。ｔは時間ステップであり、Ｎｔはシミュレーションを行う総計算時間ステップ数である。また、ｔｘは基本グリッドに基づく第１演算を行うか、サブグリッドに基づく第２演算を行うかに係る判定（グリッド変換の要否の判定）を行う時間ステップ（時間帯）である。時間ステップがｔｘでないときには、グリッドタイプの変換とそれに伴う第１演算と第２演算の切り換えは行わない。また、ｂはそれぞれのブロックに割り振られたブロック番号であり、Ｎｂは解析対象領域における総ブロック数である。
【００６８】
図６において、ステップＳ１００で、処理が開始されると、続くステップＳ１０１では、ｔ＝ｔxであるか否かが判定される。ここで、ｔ＝ｔxでないと判定されると、ステップＳ１０７に進み、シミュレーションのための伝搬計算を実行する。一方、ステップＳ１０１における判定がＹＥＳであれば、ステップＳ１０２に進み、ブロックｂの現時点でのグリッドタイプが、サブグリッドであるか否かが判定される。
【００６９】
当該ブロックｂのグリッドタイプがサブグリッドであり、ステップＳ１０２における判定がＹＥＳである場合には、ステップＳ１０３に進み、当該ブロックｂ内の場の物理量が閾値より小さいか否かが判定される。ステップＳ１０３における判定がＮＯであれば、サブグリッドに基づく第２演算から、基本グリッドに基づく第１演算に切り換えることなく、そのままステップＳ１０６に進む。
【００７０】
一方、ステップＳ１０３における判定がＹＥＳであれば、ステップＳ１０４において、サブグリッドに基づく第２演算から、基本グリッドに基づく第１演算に切り換えを行う。次の、ステップＳ１０５では、必要となる第１格子点の値を抽出し、その他の第２格子点の物理量は破棄してメモリを解放し、ステップＳ１０６に進む。
【００７１】
ステップＳ１０２において、当該ブロックｂのグリッドタイプが基本グリッドであり、判定がＮＯである場合には、ステップＳ１１０に進み、当該ブロックｂ内の場の物理量が閾値以上であるか否かが判定される。ステップＳ１１０における判定がＮＯであれば、基本グリッドに基づく第１演算から、サブグリッドに基づく第２演算に切り換えることなく、そのままステップＳ１０６に進む。
【００７２】
一方、ステップＳ１１０における判定がＹＥＳであれば、ステップＳ１１１において、基本グリッドに基づく第１演算からサブグリッドに基づく第２演算に切り換えを行う。次の、ステップＳ１１２では、第１格子点に基づいて、必要となる第２格子点の物理量を、前記したように補間により生成する。
【００７３】
ステップＳ１０６では、ｂ＝Ｎｂであるか否かが判定され、当該判定がＮＯであれば、
ステップＳ１１３に進み、ｂを１インクリメントして、再びステップＳ１０２に戻り、ループする。
【００７４】
ステップＳ１０７では、当該ブロックｂ内の各格子点について、サブグリッドＣＩＰ法による伝搬計算を実行する。
【００７５】
ステップＳ１０８では、時間ステップ更新し、ステップＳ１０９では、ｔ＞Ｎｔであるか否かを判定する。ステップＳ１０９の判定結果がＮＯである場合には、再びステップＳ１０１に戻り、ループする。
【００７６】
テップＳ１０９の判定結果がＹＥＳである場合には、全てのシミュレーション対象の総時間が経過したこととなるので、ステップＳ１１４に進み、処理を終了する。
【００７７】
以上のような、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、物理量が処理の閾値以上となる主要な波面の付近のみサブグリッドに基づく第２演算を実行することで、全体の格子点数を抑えてシミュレーションに必要なメモリ及び計算時間を大幅に低減できる。
【００７８】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、波面の移動に合わせてサブグリッドに基づく演算を行う領域を設定するので、領域全体をサブグリッドと同じ離散化幅で細かく離散化した場合と同等程度の計算精度が得られる。
【００７９】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、解析対象領域を等サイズの矩形ブロックに分割することで、ブロック境界において隣接するブロックが１体１で対応するので、サブグリッドに基づく演算への切り換え時及び境界部分の伝搬計算に必要な隣接ブロックのデータを参照するアルゴリズムが簡単になる。
【００８０】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、ブロックサイズを適切に設定することにより、様々な計算条件に柔軟に対応できる。例えば、計算時間はかかっても使用メモリを節約したい場合は、小さいブロックサイズを設定する。サブグリッド領域を最小限に設定して必要メモリを少なくすることができる。この場合、ブロック数が多いためグリッド変換やブロック境界の伝搬計算にかかる計算時間は増える。逆に、使用メモリより計算時間を節約したい場合は、大きいブロックサイズを設定する。波面付近のやや広い範囲にサブグリッド領域を設定するため必要なメモリは増えるが、ブロック数が少ないためグリッド変換やブロック境界の伝搬計算にかかる計算時間は減る。
【００８１】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムによれば、基本グリッドにより演算を行うか、サブグリッドにより演算を行うか、の判定と、その切り換えを行う時間ステップ間隔を大きく取れば、計算全体で判定と変換の回数を減らせ、グリッド変換処理による計算時間の増加を少なくすることができる。
【００８２】
また、本発明に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムは、音波の伝搬解析以外にも、電磁波など他の波動一般の伝搬解析にも適用可能である。
（実施例）
以下、２次元自由音場の解析例を示す。基本グリッドとサブグリッドの空間離散化幅は
、それぞれ０．０９ｍ、０．０３ｍとし、時間離散化幅はΔｔ＝０．４×１０-4ｓとした。初期条件は、音場の中心に振幅１及び分散０．２のＧａｕｓｓｉａｎを与えた。音場を１５×１５のブロックに分割し、各ブロックに１６×１６の第１格子点を定義した。音圧０．０２以上をサブグリッドの生成条件とした。また、比較として空間全体を０．０９ｍ或いは０．０３ｍで離散化した計算も行った。
【００８３】
図７に音圧の空間分布を、図８に音場の中心から５．７６ｍ点における音圧波形を示す。図７において、（Ａ）は解析対象領域をサブグリッドとしてシミュレーションを行った場合、（Ｂ）は本発明に係る方法でシミュレーションを行った場合を示している。また、図７（Ｂ）中における矩形（□）は、サブグリッドに基づく第２演算が行われたブロックを示している。図７（Ｂ）に示す本発明による計算結果は、解析対象領域全体を詳細に離散化した図７（Ａ）に示す結果と良く一致していることが分かる。
【００８４】
また、図８における囲み枠は、波形のピーク近傍を拡大するものであるが、本発明に係る方法によれば、解析対象領域全体を詳細に離散化（サブグリッド化）して演算を行ったものと遜色がなく、音場を粗く離散化（基本グリッド化）した場合に生じる数値散逸を回避できていることが確認できる。
【００８５】
図９には、時間ステップ毎の総グリッド数を示す。時間と共に波面が広がるため、サブグリッド領域の増加に伴い総グリッド数が増えている。しかしながら、最大でも音場全体（解析対象領域全体）を詳細に離散化（サブグリッド化）した場合の１／２以下となっている。また、本発明の方法と比較例に係る方法の計算時間を表１に示す。
【００８６】
【表１】
2013214299.tif 000014
計算環境にはスクリプト言語であるＭＡＴＬＡＢを使用したため暫定的な結果であるが、本発明による方法を用いた場合の計算時間は音場全体（解析対象領域全体）を詳細に離散化（サブグリッド化）した場合に比べ約１／５である。
【００８７】
次に、本発明の第２実施形態について説明する。本発明の基本的な考え方は、解析対象領域における物理量をシミュレーションする際に、解析対象領域全体を複数のブロックに分割し、それぞれのブロック内の物理量の値に応じてブロック内を細かく離散化して解析する場合（サブグリッドによる第２演算を行う場合）と粗く離散化して解析する場合（基本グリッドによる第１演算を行う場合）を計算過程で切り換えることにより、使用メモリを減らし計算時間を短くする、というものである。
【００８８】
第２実施形態においては、さらに、基本グリッドにより演算を行うか、サブグリッドにより演算を行うか、の判定と、その切り換えを行う時間ステップ間隔を大きく取ることで、計算全体で判定と変換の回数を減し、グリッド変換処理による計算時間の増加を少なくすることを目的としている。
【００８９】
まず、これまでに説明した本発明に係る物理量シミュレーション方法の考え方の概略について、再度、図面を参照して説明する。図１０は本発明に係る物理量シミュレーション方法の概念を示す図である。
【００９０】
図１０は、中心０で発生した音の音圧波形の伝搬をシミュレーションしている様子を示している。図１０において、音圧波形が伝搬する先頭付近のブロック（実線で示すブロック）では、ブロック内の場の物理量は閾値以上となり、サブグリッドによる第２演算を行う。逆に、音圧波形が到達していないブロック（点線で示すブロック）、及び、音圧波形の先頭付近が伝搬し終わった後のブロック（点線で示すブロック）では、ブロック内の場の物理量は閾値より小さくなるので、基本グリッドによる第２演算を行う。
【００９１】
これに対して、第２実施形態においては、ブロックｂ内の場の物理量≧閾値となり、サブグリッドにより第２演算を行うように切り換えるブロックに加え、それに隣接するブロックについても、サブグリッドによる第２演算を適用する。図１１は第２実施形態に係る物理量シミュレーション方法の概念を示す図である。
【００９２】
図１１は、図１０と同じ前提を図示したものであるが、図１１においては、第１の実施形態でサブグリッドによる第２演算を行っていたブロック（実線で示すブロック）に、隣接するブロック（実線に×印が付されたブロック）についても、サブグリッドによる第２演算を行うことを示している。
【００９３】
第１実施形態に比べて、他の実施形態においては、サブグリッドによる第２演算を行うブロックは増えることとなるが、基本グリッドにより演算を行うか、サブグリッドにより演算を行うか、の判定と、その切り換えを行う時間ステップ間隔を大きく取ることで、計算全体で判定と変換の回数を減し、グリッド変換処理による計算時間の増加を少なくするようにしている。
【００９４】
音圧波形の先頭付近である判定対象の物理量が閾値を越える領域は計算の進展に伴い移動していくが、時間Ｔb＝Ｌb／ｃの範囲内であれば、その領域は、第１の実施形態でサブグリッドによる第２演算を行うブロック（実線で示すブロック）、及び、それに隣接するブロック（実線に×印が付されたブロック）内に確実に留まっている。ここで、Ｌbはブ
ロック１辺の長さであり、ｃは音速である。
【００９５】
なお、音速ｃを一般化すると、解析対象領域における媒質中で波動が伝搬する速さ、とすることができる。したがい、例えば、本実施形態を、音波の伝搬解析以外の、電磁波などの波動の伝搬解析に適用すると、ｃとしては光速を用いるようにする。
【００９６】
以上のような、第２実施形態においては、一度基本グリッドによる演算を行うかサブグリッドにより演算を行うかの判定と演算の切り換えを行った後は、時間Ｔbを経過するま
では再度の判定は必要なく、判定と切り換えを行う時間ステップ間隔を大きくとることができる。
【００９７】
ここで、第２実施形態に係る物理量シミュレーション方法をコンピューターに実行させる際のアルゴリズムについて説明する。図１２は第２実施形態に係る物理量シミュレーション方法をコンピューターに実行させるためのフローチャート例を示す図である。
【００９８】
以下、図７に示した第１実施形態のフローチャートと相違する点について、説明する。
【００９９】
第２実施形態に係る物理量シミュレーション方法では、時間Ｔb経過毎に判定を行うこ
とが前提であるので、ステップＳ１０１と異なり、ステップＳ２０１ではｔ＝ｎＴb（た
だし、nは自然数）であるか否かが判定され、ステップＳ２０１における判定がＹＥＳで
あるときのみ、ステップＳ２０２に進むこととなる。
【０１００】
また、第２実施形態に係る物理量シミュレーション方法では、第１の実施形態でサブグリッドによる第２演算を行うブロック（実線で示すブロック）に加え、それに隣接するブロック（実線に×印が付されたブロック）についても、演算方法を切り換えるので、ステップＳ１１１と異なり、ステップＳ２１１で、ブロックｂと、ブロックｂに隣接するブロックをサブグリッドに基づく第２演算に切り換える。
【０１０１】
以上のような、第２実施形態に係る物理量シミュレーション方法によれば、基本グリッドにより演算を行うか、サブグリッドにより演算を行うか、の判定と、その切り換えを行う時間ステップ間隔を大きく取ることで、計算全体で判定と変換の回数を減し、グリッド変換処理による計算時間の増加を少なくすること可能となる。
【０１０２】
なお、サブグリッドに基づく第２演算に切り換える範囲を、上記当該ブロック、その隣接ブロック、更に隣接するブロックに隣接するブロックにまで拡げれば、一度判定と切り換えを行った後、時間２Ｔbを経過するまでは再度の判定は必要なく、判定と切り換えを
行う時間ステップ間隔を更に大きくとることができる。以下同様に、第２演算に切り換える範囲を拡げるほど、判定と切り換えを行う時間ステップ間隔は大きくとることができる。
【０１０３】
以上、第２実施形態によれば、可変グリッドＣＩＰ法による波動伝搬シミュレーションにおいて、基本グリッドによる演算を行うかサブグリッドにより演算を行うか、の判定と、その切り換えを行う時間ステップ間隔を大きくとることができる。また、前記判定とグリッド変換処理の回数を減らすことにより、全体の計算時間の増加を少なくすることができる。
【０１０４】
また、第２実施形態に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムは、音波の伝搬解析以外にも、電磁波など他の波動一般の伝搬解析にも適用可能である。
【０１０５】
次に、本発明の第３実施形態について説明する。図１３は第３実施形態に係る物理量シミュレーション方法の概念を示す図である。
【０１０６】
図１１に示した第２実施形態では、ブロックｂ内の場の物理量≧閾値となり、サブグリッドにより第２演算を行うように切り換えるブロックに加え、それに隣接するブロックについても、サブグリッドによる第２演算を適用するようにしていた。しかしながら、これによれば、ブロック内の場の物理量が閾値より小さくなっている可能性が高い、音圧波形の先頭付近が伝搬し終わった後のブロックについても、サブグリッドによる第２演算を適用することとなる。
【０１０７】
そこで、これを避けるために、第３実施形態においては、ＣＩＰ法により物理量を演算する際に利用される、特定方向に伝搬する波動の成分を表す特性値ｆ±、ｇ±を利用した判定を行うようにしている。
【０１０８】
図１３は、図１１と同じ前提を図示したものであるが、図１１においては、音圧波形の先頭付近が伝搬し終わった後のブロックもサブグリッドによる第２演算を適用されているが、図１３に示す第３実施形態では、そのようなことがない。
【０１０９】
第３実施形態は、判定対象となる物理量が音圧、各軸方向の粒子速度、それらの方向微分値から算出した、ｆ±、ｇ±（特定方向に伝搬する波動の成分を表す特性値）のいずれかの場合に適用される。
【０１１０】
ｆ+は正方向、ｆ-は負方向に伝搬する波動の成分を表し、ｇ±はそれらの方向微分値である。なお、以下では２次元解析を例とし、±ｘ方向に伝搬する成分をｆx±、±ｙ方向
に伝搬する成分をｆy±、それらの方向微分値をｇx±、ｇy±として表す。
【０１１１】
具体的には、あるブロック内で判定対象の物理量（ｆx±、ｆy±、又はｇx±、ｇy±）が設定した閾値を越える場合、当該ブロックとその物理量の伝搬方向にあたる隣接ブロックをサブグリッドに基づく第２演算に切り換える。
【０１１２】
例えば、ｆx±、ｆy±を判定対象とした場合、ｆx+が閾値を越えていれば当該ブロックと＋ｘ方向に隣接するブロックを第２演算に切り換える。或いは、ｆy―が閾値を越えて
いれば当該ブロックと−ｙ方向に隣接するブロックを第２演算に切り換える。
【０１１３】
また、ｆx+およびｆy―が閾値を越えていれば、＋ｘに隣接するブロック、かつ、−ｙ
に隣接するブロック、すなわち、当該ブロックの斜め方向に隣接するブロックについても第２演算に切り換える。
【０１１４】
音圧波形の先頭付近である判定対象の物理量が閾値を越える領域は計算の進展に伴い移動していくが、時間Ｔb＝Ｌb／ｃの範囲内であれば、その領域は上記のような方法によりサブグリッドに変換した上記当該ブロック或いは隣接ブロック内に確実に留まっている。ここで、Ｌbはブロック１辺の長さであり、ｃは音速である。
【０１１５】
以上のような、第３実施形態においては、一度基本グリッドによる演算ステップを行うかサブグリッドにより演算を行うかの判定と演算の切り換えを行った後は、時間Ｔbが経
過するまでは再度の判定は必要なく、判定と切り換えを行う時間ステップ間隔を大きくとることができる。
【０１１６】
第３実施形態によれば、波動の伝搬方向を考慮することにより、第２実施形態と比較して第２演算に切り換える範囲を限定することで使用メモリと計算時間を減らすことができる。
【０１１７】
ここで、第３実施形態に係る物理量シミュレーション方法をコンピューターに実行させる際のアルゴリズムについて説明する。図１４は第３実施形態に係る物理量シミュレーション方法をコンピューターに実行させるためのフローチャート例を示す図である。
【０１１８】
以下、図７に示した第１実施形態のフローチャートと相違する点について、説明する。
【０１１９】
第３実施形態に係る物理量シミュレーション方法では、時間Ｔb経過毎に判定を行うこ
とが前提であるので、ステップＳ１０１と異なり、ステップＳ３０１ではｔ＝ｎＴb（た
だし、nは自然数）であるか否かが判定され、ステップＳ３０１における判定がＹＥＳで
あるときのみ、ステップＳ３０２に進むこととなる。
【０１２０】
第３実施形態に係る物理量シミュレーション方法では、判定対象が、特定方向に伝搬する波動の成分を表す特性値（ｆx±、ｆy±、又はｇx±、ｇy±）であるので、ステップＳ１０３と異なり、ステップＳ３０３では、「特性値＜閾値」が判定される。同様に、ステップＳ１１０と異なり、ステップＳ３１０では、「特性値≧閾値」が判定される。
【０１２１】
また、第３実施形態に係る物理量シミュレーション方法では、第１の実施形態でサブグリッドによる第２演算を行うブロック（実線で示すブロック）に加え、特定方向でそれに隣接するブロック（実線に×印が付されたブロック）についても、演算方法を切り換えるので、ステップＳ１１１と異なり、ステップＳ３１１で、ブロックｂと、特定方向でブロ
ックｂに隣接するブロックをサブグリッドに基づく第２演算に切り換える。
【０１２２】
以上のような、第３実施形態に係る物理量シミュレーション方法によれば、基本グリッドにより演算を行うか、サブグリッドにより演算を行うか、の判定と、その切り換えを行う時間ステップ間隔を大きく取ることで、計算全体で判定と変換の回数を減し、グリッド変換処理による計算時間の増加を少なくすること可能となる。
【０１２３】
なお、サブグリッドに基づく第２演算に切り換える範囲を、上記当該ブロック、その隣接ブロック、更に隣接ブロックに隣接するブロックにまで拡げれば、一度判定と切り換えを行った後、時間２Ｔbを経過するまでは再度の判定は必要なく、判定と切り換えを行う
時間ステップ間隔を更に大きくとることができる。以下同様に、第２演算に切り換える範囲を拡げるほど、判定と切り換えを行う時間ステップ間隔は大きくとることができる。
【０１２４】
以上、第３実施形態によれば、可変グリッドＣＩＰ法による波動伝搬シミュレーションにおいて、基本グリッドによる演算を行うかサブグリッドにより演算を行うか、の判定と、その切り換えを行う時間ステップ間隔を大きくとることができる。また、前記判定とグリッド変換処理の回数を減らすことにより、全体の計算時間の増加を少なくすることができる。
【０１２５】
また、第３実施形態に係る物理量シミュレーション方法及びそれを用いた物理量シミュレーションシステムは、音波の伝搬解析以外にも、電磁波など他の波動一般の伝搬解析にも適用可能である。
【符号の説明】
【０１２６】
１０・・・システムバス
１１・・・ＣＰＵ（Ｃｅｎｔｒａｌ  Ｐｒｏｃｅｓｓｉｎｇ  Ｕｎｉｔ）
１２・・・ＲＡＭ（Ｒａｎｄｏｍ  Ａｃｃｅｓｓ  Ｍｅｍｏｒｙ）
１３・・・ＲＯＭ（Ｒｅａｄ  Ｏｎｌｙ  Ｍｅｍｏｒｙ）
１４・・・通信制御部
１５・・・入力制御部
１６・・・出力制御部
１７・・・外部記憶装置制御部
１８・・・入力部
１９・・・出力部
２０・・・外部記憶装置
２１・・・インターフェイス部
２１・・・グラフィック制御部
２２・・・ディスプレイ装置
【図１】
2013214299.tif 000015
【図２】
2013214299.tif 000016
【図３】
2013214299.tif 000017
【図４】
2013214299.tif 000018
【図５】
2013214299.tif 000019
【図６】
2013214299.tif 000020
【図７】
2013214299.tif 000021
【図８】
2013214299.tif 000022
【図９】
2013214299.tif 000023
【図１０】
2013214299.tif 000024
【図１１】
2013214299.tif 000025
【図１２】
2013214299.tif 000026
【図１３】
2013214299.tif 000027
【図１４】
2013214299.tif 000028