JP2023170345A - Ｆｄｔｄ解析方法 - Google Patents

Abstract

【課題】 電磁波および音波の伝搬を解析するＦＤＴＤ法において精度向上および計算時間短縮の両立を図る。【解決手段】 電磁波および音波の伝搬を解析するＦＤＴＤ法において、奇数階微分項を付加した差分式を用いる。６次精度のＦＤＴＤ法の場合は、４次精度の３階微分項×α、２次精度の５次微分項×βを付加した差分式とすることができる。係数α、βは、クーラン数Ｃと相関を有しており、分散関係式によって伝搬する波の位相速度を求め、理論値との誤差が最小となるように探索的に設定するこうして得られた差分式を用いて解析することにより、精度向上を図りつつ、クーラン数を大きくすることができ、演算の所要時間を短縮することができる。【選択図】 図３

Description

本発明は、ＦＤＴＤ（Finite-Difference Time-Domain）法による解析方法に関する。

電磁波伝搬及び音波伝搬の数値解析の方法として、ＦＤＴＤ法が広く用いられている。ＦＤＴＤ法について電磁波を例にとって概要を説明する。

電磁波の伝搬を記述するＭａｘｗｅｌｌ方程式は次式で与えられる。

ＦＤＴＤ法は、Ｍａｘｗｅｌｌ方程式を時間及び空間差分の形で表したものであり、第一式、即ち磁界Ｂについて２次精度で近似した式は次式となる。

また第二式、即ち電解Ｅについて２次精度で表した式は次式となる。

ここで、Ｅ、Ｂの右上に示したｔ等の符号は、時刻ｔにおける電界（Ｅ）、磁界（Ｂ）であることを意味し、Ｅ、Ｂの右下に示したｘ、ｙ、ｚは、それぞれ電界強度及び磁界強度の方向を表している。また、Δｘ、Δｙ、Δｚは、それぞれｘ、ｙ、ｚ方向の空間格子の刻み幅を表し、Δｔは時間積分の刻み幅を表している。
上記２つの式において、ｘ、ｙ、ｚをそれぞれｙ、ｚ、ｘと置き換えることによりｙ成分についての関係式を得ることができる。また、ｘ、ｙ、ｚをそれぞれｚ、ｘ、ｙと置き換えることによりｚ成分についての関係式を得ることができる。
２次精度とは、Δｔ、Δｘ、Δｙ、Δｚの２乗の精度で表した式であることを意味する。

Ｍａｘｗｅｌｌ方程式の第一式を４次精度で表した次式も知られている。

同じく第二式を４次精度で表した式は次式となる。

ＦＤＴＤ法については、種々の改良が提案されており、例えば、特許文献１は、解析対象となる領域を複数のブロックに分割して演算することにより、解析の所要時間を短縮する方法を提案している。

K. S. Yee, Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic media, IEEE Transactions on Antennas and Propagation, AP-14 (1966) 302-307. P. G. Petropoulos, Phase error control for FD-TD methods of second and fourth order accuracy, IEEE Transactions on Antennas and Propagation, 42 (1994) 859-862. E. Ikata and G. Tay, Finite-difference time domain acoustic-wave algorithm, Il Nuovo Cimento D, 20 (1998) 1779-1793.

特開２００６－１３９７２３号公報

ＦＤＴＤ法では、数値不安定を生じることなく演算を可能とするための時間刻みΔtあたりの電磁波あるいは音波の伝搬速度の上限値が決まっており、これをクーラン条件とよぶ。空間格子の刻み幅が等しい場合（Δｘ＝Δｙ＝Δｚ）を仮定すると、２次精度のＦＤＴＤ法の伝搬速度の上限値を表すクーラン数Ｃは、３次元では以下の式で与えられる。

また、４次精度の場合は、クーラン数は次式で与えられる。

このように、空間差分の精度を上げると、時間刻みの上限値は小さくなるため、電磁波あるいは音波の伝搬を所定の時間まで解析するために必要な時間ステップ数（時間積分の回数）が増える。すなわち、空間差分の精度を上げることによって計算時間が非常に長くなるという課題があった。
本発明は、この課題に鑑み、ＦＤＴＤ法による解析において、計算精度の向上と計算時間の短縮を両立させることを目的とする。

本発明は、
コンピュータによって、ＦＤＴＤ法により電波または音波の伝搬を表す２つの方程式に含まれる２つの状態量を数値計算するＦＤＴＤ解析方法であって、
前記コンピュータが実行する工程として、
（ａ）空間内に設定された各格子点における２つの前記状態量その他の初期条件を入力し、前記コンピュータのメモリの所定領域に格納するステップと、
（ｂ）前記方程式の第一式を表す差分式によって、所定の時間刻みΔｔ後の前記格子点における第１の状態量を計算し、結果を前記コンピュータのメモリの所定領域に格納するステップと、
（ｃ）前記方程式の第二式を表す差分式によって、所定の時間刻みΔｔ後の前記格子点における第２の状態量を計算し、結果を前記コンピュータのメモリの所定領域に格納するステップと、
（ｄ）所定の終了条件に至るまで前記ステップ（ｂ）、（ｃ）を繰り返し実行するステップとを備え、
前記差分式は、ＦＤＴＤ法の計算式に、奇数階微分に係数を乗じた付加項を１つ以上加算した式であるＦＤＴＤ解析方法。

発明者は、ＦＤＴＤ解析方法について研究を重ねた結果、上述の形の式を用いることにより、係数およびクーラン数を調整すれば、精度向上と計算速度の短縮の両立を図ることができることを見いだした。上述の式は、semi-Lagrange法の考え方をＦＤＴＤ解析に応用することにより着想されたものである。ただし、semi-Lagrange法では偶数階微分項も付加するのに対し、偶数階微分項には波形を減衰させる数値散逸効果があることを考慮し、本発明では奇数階微分項のみを付加している点に特徴を有する。
ＦＤＴＤ法の精度および付加項の精度は種々の組み合わせが可能である。また、付加項は１つだけである必要はなく、２以上の項を付加してもよい。

本発明において、
前記差分式は、ＦＤＴＤ法におけるｎ（ｎは偶数）次精度の計算式に、（ｎ－２ｍ）（ｍはｎ／２より小さい自然数）次以下の精度の（２ｍ＋１）階微分に係数を乗じた付加項を１つ以上加算した式としてもよい。

研究の結果、差分式に用いる点数を考慮すると、上述した形式が最も効率が良いことが見いだされた。ｎ＝４とすれば、差分式は「４次精度のＦＤＴＤ法＋２次精度の３階微分」となる。ｎ＝６とすれば、差分式は「６次精度のＦＤＴＤ法＋４次精度の３階微分＋２次精度の５階微分」となる。他にも種々の差分式を設定可能である。
なお、本発明の差分式は、上述した形式に限定されるものではなく、例えば、６次精度ＦＤＴＤ法＋２次精度３階微分項、２次精度ＦＤＴＤ法＋４次精度３階微分項などの形としてもよい。

本発明において、
前記付加項に乗ぜられるそれぞれの係数は、クーラン数Ｃと相関を有することが好ましい。

研究の結果、係数とクーラン数とを相関を有する状態で調整することにより、精度と計算速度を両立できることが見いだされた。両者の相関は、関数で表すことは困難であり、具体的な数の組み合わせまたはテーブルの形式で与えられることになる。
本発明による解析方法を利用する場合には、これらの相関に基づいて、演算における空間の格子および時間刻みを、数値不安定が生じない範囲で設定して差分式を特定した上で演算を行うことになる。

上述の通り係数およびクーラン数に相関を有する場合、
前記ステップ（ａ）において、前記係数またはクーラン数Ｃのいずれかを入力し、
前記ステップ（ｂ）および前記ステップ（ｃ）の実行前に、前記係数とクーラン数との間の予め設定された相関を利用して、残余の前記係数およびクーラン数を設定するステップを備えるものとしてもよい。

上記態様によれば、例えば、空間の格子および時間刻みなどの条件を設定すれば、クーラン数が決まるため、これに応じて適切な係数の値が設定され、演算を行うことができる。逆に、係数を設定すれば、これに応じてクーラン数を設定して演算を行うことができる。
予め設定された相関は、関数、テーブルなど種々の形式で用意することができる。

本発明において、
前記付加項に乗ぜられるそれぞれの係数は、前記伝搬における位相速度の誤差を最小とするという評価基準に基づいて設定されていることが望ましい。

本来であれば、電磁波または音波は全波数において一定の位相速度で波形が伝搬することが理想的であるが、数値シミュレーションにおいては波数によって異なる位相速度で伝搬するという現象が現れる。上記態様によれば、この位相速度の誤差が最小となることを評価基準として係数を設定するため、計算精度を向上させることができる。
なお、「評価基準として」とは、厳密に誤差最小を実現する係数でなくても良いことを意味する。例えば、複数の係数候補について位相速度を評価した上で、その中で誤差が最小のものを選択する方法をとってもよい。

先に説明した通り、本発明においては、ＦＤＴＤ法の精度および付加項には種々の組み合わせが可能である。
例えば、
前記差分式は、ＦＤＴＤ法における６次精度の計算式に、４次精度の３階微分に係数を乗じた付加項を加算した式としてもよい。
この場合、
前記差分式は、さらに、２次精度の５階微分に係数を乗じた付加項を加算した式としてもよい。

また、前記差分式は、ＦＤＴＤ法における４次精度の計算式に、４次精度の３階微分に係数を乗じた付加項を加算した式としてもよい。
さらに、前記差分式は、ＦＤＴＤ法における４次精度の計算式に、２次精度の３階微分に係数を乗じた付加項を加算した式としてもよい。

発明者の研究成果によれば、少なくともこれらの組み合わせにおいては、適切な係数およびクーラン数が設定でき、計算精度および計算速度を両立できることが確認されている。これらの態様における具体的な差分式の例および係数、クーラン数の組み合わせについては後述する。

本発明は、上述したＦＤＴＤ解析方法としての態様に限らず、かかる解析を実現するための解析装置、コンピュータによってかかる解析を実現するためのコンピュータプログラム、およびかかるプログラムを記録したコンピュータ読み取り可能な記録媒体として構成することもできる。ここでコンピュータ読み取り可能な記録媒体としては、フラッシュメモリ、ハードディスク、光ディスクなど種々の媒体を用いることができる。

解析装置の構成を示した説明図である。 解析処理のフローチャートである。 ３次元における伝搬速度の誤差の計算結果を示す説明図である。 ２次元における伝搬速度の誤差の計算結果を示す説明図である。 ２次元における電磁波の伝搬の様子を示す説明図である。

Ａ．解析装置の構成：
図１は、解析装置の構成を示した説明図である。スーパーコンピュータに、図示する各機能を実現するコンピュータプログラムをインストールすることによって構成した。解析装置は、複数のコンピュータ、サーバ等をネットワークで接続して構成してもよい。
以下、電磁波の伝搬を解析する場合を例にとって説明する。

以下、図中の各機能ブロックについて説明する。
入力部１１は、運動解析に要する種々の条件を入力する。例えば、空間格子の設定および各格子点における磁界および電界の初期条件、時間刻みなどの条件が挙げられる。

算出結果保持部１２は、入力部１１から入力された種々の初期条件を保持するとともに、電磁波の伝搬の解析を行う過程で算出される種々の計算結果を一時的に保持するバッファである。
磁界演算部１４は、与えられた条件下で後述するＦＤＴＤ法の差分式を用いて磁界の変化を演算する。計算結果は、逐次、算出結果保持部１２に格納しながら計算を進める。
電界演算部１５は、与えられた条件下で後述するＦＤＴＤ法の差分式を用いて電界の変化を演算する。計算結果は、逐次、算出結果保持部１２に格納しながら計算を進める。
磁界演算部１４、電界演算部１５は、別の機能部として示したが、両者を統合して構成しても差し支えない。
出力部１３は、以上で得られた解析結果を出力する。出力形式は問わない。

係数等設定部１６は、磁界演算部１４および電界演算部１５の演算に用いられる係数を設定する。後述する通り、本実施例では、磁界および電界の演算には、ＦＤＴＤ法の差分式に、所定の微分項に係数を乗じた付加項を加算した差分式を用いている。係数等設定部１６は、この付加項に乗ぜられる係数を設定するのである。
係数は、後述する通り、ＦＤＴＤ法を演算する際のクーラン数と相関を有している。クーラン数Ｃは、ｃ（光速）と時間刻みΔｔ、空間の格子間隔Δｘ（＝Δｙ＝Δｚとする）によってＣ＝ｃ×Δｔ／Δｘで算出される。即ち、空間の格子および時間刻みが初期条件として指定されればクーラン数は定まることになる。そして、クーラン数が定まれば、それと相関のある係数も定まることになる。
本実施例では、種々のクーラン数に対して係数を求めた結果を設定テーブル１７として予め用意しておいた。従って、係数等設定部１６は、クーラン数に応じて、設定テーブル１７を参照することにより係数を定めることができるのである。
なお、差分式で用いる係数等を初期条件の一つとしてオペレータが入力する方法をとってもよい。この場合は、係数等設定部１６および設定テーブル１７は省略することができる。

Ｂ．解析処理：
図２は、解析処理のフローチャートである。図１に示した解析装置が実行する処理である。
解析処理を開始すると、解析装置は、初期条件を入力する（ステップＳ１）。入力される解析条件は、入力部１１（図１参照）において説明した通りである。
次に、解析装置は、差分式の係数等を算出する（ステップＳ２）。係数等の算出方法は、係数等設定部１６（図１参照）において説明した通りである。係数も入力する場合には、ステップＳ２の処理は省略してもよい。

そして、解析装置は、差分式を用いて磁界を演算し、その結果を格納する（ステップＳ３）。また、差分式を用いて電界を演算し、その結果を格納する（ステップＳ４）。上記演算を終了条件が満たされるまで（ステップＳ５）、繰り返し実行して、結果を出力し（ステップＳ６）、解析処理を終了する。
終了条件としては、予め設定した時間後の結果が得られた場合、計算時間が所定の制限値を超えた場合、その他エラーが生じた場合などとすることができる。

Ｃ．差分式：
以下、実施例の解析において用いられる差分式について具体例をあげて説明する。本実施例では、ＦＤＴＤ法における計算式に、奇数階微分に係数を乗じた付加項を１つ以上加算した式を用いている。付加項は、１つでも２以上でも良い。
このように差分式は、精度および付加項の選択によって、多様な組み合わせが可能であるが、以下では、２次元および３次元の場合に、それぞれ係数とクーラン数の相関の検討結果が得られているものを具体例として示す。

Ｃ１．３次元６次精度＋４次精度３階微分項＋２次精度５階微分項：
かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第一式、即ち磁界を表す差分式は、次の通りとなる。

ここで、Ｔ１（４，３）は、４次精度３階微分項を表し、次式で与えられる。次式におけるαは係数である。

また、Ｔ２（２，５）は、２次精度５階微分項を表し、次式で与えられる。次式におけるβは係数である。

また、Ｍａｘｗｅｌｌ方程式の第二式、即ち電界を表す差分式は、次の通りとなる。

ここで、Ｔ３（４，３）は、４次精度３階微分項を表し、次式で与えられる。次式におけるαは係数である。

また、Ｔ４（２，５）は、２次精度５階微分項を表し、次式で与えられる。次式におけるβは係数である。

以上の式では、ｘ成分のみを示しているが、各式において、ｘ、ｙ、ｚをそれぞれｙ、ｚ、ｘと置き換えることによりｙ成分についての関係式を得ることができる。また、ｘ、ｙ、ｚをそれぞれｚ、ｘ、ｙと置き換えることによりｚ成分についての関係式を得ることができる。３次元については、以下同様である。

付加項の係数α、βは、まず数値不安定を生じないよう設定する必要がある。ＦＤＴＤ法の差分式をフーリエ変換することで、周波数ωとクーラン数、時間刻Δｔ、格子間隔Δｘなどの関数として波動の分散関係式が得られるが、この分散関係式において周波数ωが虚数解となるとき、数値不安定が生じることになる。従って、係数α、βは、分散関係式において周波数ωが虚数解とならないことを条件の一つとして設定する必要がある。
また、電磁波は全波数において一定の位相速度で波形が伝搬することが理想的であるため、分散関係式から得られる位相速度、即ちωの実部を波数κで除した値、と理論値の２乗平均誤差をκについて積分した値が最も小さくなるよう係数を設定することが好ましい。本実施例では、これらの条件に基づいて、クーラン数および係数の組み合わせを探索的に求めることとした。以下で説明する種々の差分式においても同様である。

本ケースの場合、分散関係式は次式で与えられる。

空間格子がｘ、ｙ、ｚ方向に等しい間隔（Δｘ＝Δｙ＝Δｚ）の場合、この分散関係式に基づいて得られるクーラン数および係数の組み合わせは、次の表の通りである。Δｘ、Δｙ、Δｚが等しくない場合は、ｘ方向、ｙ方向、ｚ方向のそれぞれに対して、下表を求めることになる。

図３は３次元における伝搬速度の誤差の計算結果を示す説明図である。横軸にクーラン数をとり、縦軸に伝搬速度についての演算結果と理論値との誤差を示した。ＦＤＴＤ（２，２）は従来の２次精度のＦＤＴＤ法による演算結果を示している。クーラン条件を満たす範囲、即ちクーラン数Ｃが約０．５７７より小さい範囲でしか解が得られない。ＦＤＴＤ（２，４）は、従来の４次精度のＦＤＴＤ法による演算結果を示している。同じくクーラン条件を満たす範囲、即ちクーラン数Ｃが約０．４９４より小さい範囲でしか解が得られない。

Ｐｒｅｓｅｎｔは、本実施例における結果である。従来のＦＤＴＤ法による結果とは異なり、クーラン数Ｃが１．０以下の範囲で解が得られている。クーラン数が１．０より大きい範囲であっても、それに応じた適切なα、βを探索できれば解を得ることは可能である。このように、本実施例の方法によれば、従来のＦＤＴＤ法よりも、クーラン数を大きくすることができる。クーラン数は、ｃΔｔ／Δｘで求まるから、空間格子の間隔Δｘが一定であれば時間刻みΔｔを大きくすることができることを意味する。
本実施例の方法では、付加項の分、１回の演算量が増大するものの、本実施例によれば、従来のＦＤＴＤ法よりも時間刻みを大きくすることができ、全体として計算時間を従来の２次精度のＦＤＴＤ法と同等にすることが可能となる。

また、従来のＦＤＴＤ法で解が得られている範囲で比較すると、本実施例の演算結果は、従来のＦＤＴＤ法よりも全体に誤差が小さくなっていることがわかる。即ち、本実施例の方法によれば、従来のＦＤＴＤ法よりも精度を向上させることが可能となる。

Ｃ２．２次元６次精度＋４次精度３階微分項＋２次精度５階微分項：
かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第一式、即ち磁界を表す差分式は、ｘ成分について次の通りとなる。

また、ｙ成分について次の通りとなる。

かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第二式、即ち電界を表す差分式は、次の通りとなる。

以上の３式において、係数αが乗ぜられた項が４次精度３階微分の付加項であり、係数βが乗ぜられた項が２次精度５階微分の付加項である。
本ケースの場合、係数α、βおよびクーラン数を探索するための分散関係式は次式で与えられる。

空間格子がｘ、ｙ方向に等しい間隔（Δｘ＝Δｙ）の場合、この分散関係式に基づいて得られるクーラン数および係数の組み合わせは、次の表の通りである。Δｘ、Δｙが等しくない場合は、ｘ方向、ｙ方向のそれぞれに対して、下表を求めることになる。

図４は２次元における伝搬速度の誤差の計算結果を示す説明図である。横軸にクーラン数をとり、縦軸に伝搬速度についての演算結果と理論値との誤差を示した。ＦＤＴＤ（２，２）は従来の２次精度のＦＤＴＤ法による演算結果、ＦＤＴＤ（２，４）は従来の４次精度のＦＤＴＤ法による演算結果、Ｐｒｅｓｅｎｔは本実施例による演算結果を示している。３次元の場合と同様、本実施例によれば、クーラン数を大きくすることができるため、
全体として計算時間を短縮することが可能となることがわかる。また、従来のＦＤＴＤ法よりも精度を向上させることも可能となる。

図５は、２次元における電磁波の伝搬の様子を示す説明図である。初期条件として原点（ｘ＝０，ｙ＝０）でピークを有するガウス関数を与えた場合の所定時間経過後の電磁波の伝搬の結果を表している。
図５（ａ）は、従来の２次精度のＦＤＴＤ法をクーラン数０．５で実行した結果である。図５（ｂ）は、従来の４次精度のＦＤＴＤ法をクーラン数０．５で実行した結果である。図５（ｃ）は、本実施例の方法をクーラン数０．５で実行した結果である。図５（ｄ）は本実施例の方法を、クーラン数１．０で実行した結果である。
図示する通り、クーラン数０．５での本実施例の結果（図５（ｃ））は、従来の２次精度（図５（ａ））および４次精度（図５（ｂ））とほとんど差違のない結果が得られている。クーラン数１．０での本実施例の結果（図５（ｄ））では、ｘ方向およびｙ方向に、若干の数値振動による縞模様が現れているが、それ以外は、２次精度（図５（ａ））および４次精度（図５（ｂ））と大きな差違はない。即ち、本実施例の方法によれば、クーラン数を１．０まで大きくしても、看過できないほどの誤差を与えることなく演算結果が得られることがわかる。

Ｃ３．３次元６次精度＋４次精度３階微分項：
かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第一式、即ち磁界を表す差分式は、次の通りとなる。

また、Ｍａｘｗｅｌｌ方程式の第二式、即ち電界を表す差分式は、次の通りとなる。

以上の２式において、係数αが乗ぜられた項が４次精度３階微分の付加項である。
本ケースの場合、係数αおよびクーラン数を探索するための分散関係式は次式で与えられる。

空間格子がｘ、ｙ、ｚ方向に等しい間隔（Δｘ＝Δｙ＝Δｚ）の場合、この分散関係式に基づいて得られるクーラン数および係数αの組み合わせは、次の表の通りである。Δｘ、Δｙ、Δｚが等しくない場合は、ｘ、ｙ、ｚ方向のそれぞれに対して、下表を求めることになる。

既に説明したケースと同様、本ケースにおいても、計算精度を確保または向上しつつ、クーラン数を大きく設定することにより演算の所要時間を短くすることができる。

Ｃ４．２次元６次精度＋４次精度３階微分項：
かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第一式、即ち磁界を表す差分式は、ｘ成分について次の通りとなる。

また、ｙ成分について次の通りとなる。

かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第二式、即ち電界を表す差分式は、次の通りとなる。

以上の各式において、係数αが乗ぜられた項が４次精度３階微分の付加項である。
本ケースの場合、係数αおよびクーラン数を探索するための分散関係式は次式で与えられる。

空間格子がｘ、ｙ方向に等しい間隔（Δｘ＝Δｙ）の場合、この分散関係式に基づいて得られるクーラン数および係数の組み合わせは、次の表の通りである。Δｘ、Δｙが等しくない場合は、ｘ方向、ｙ方向のそれぞれに対して、下表を求めることになる。

既に説明したケースと同様、本ケースにおいても、計算精度を確保または向上しつつ、クーラン数を大きく設定することにより演算の所要時間を短くすることができる。

Ｃ５．３次元４次精度＋４次精度３階微分項：
かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第一式、即ち磁界を表す差分式は、次の通りとなる。

また、Ｍａｘｗｅｌｌ方程式の第二式、即ち電界を表す差分式は、次の通りとなる。

以上の２式において、係数αが乗ぜられた項が４次精度３階微分の付加項である。
本ケースの場合、係数αおよびクーラン数を探索するための分散関係式は次式で与えられる。

空間格子がｘ、ｙ、ｚ方向に等しい間隔（Δｘ＝Δｙ＝Δｚ）の場合、この分散関係式に基づいて得られるクーラン数および係数αの組み合わせは、次の表の通りである。Δｘ、Δｙ、Δｚが等しくない場合は、ｘ方向、ｙ方向、ｚ方向のそれぞれに対して、下表を求めることになる。

既に説明したケースと同様、本ケースにおいても、計算精度を確保または向上しつつ、クーラン数を大きく設定することにより演算の所要時間を短くすることができる。

Ｃ６．２次元４次精度＋４次精度３階微分項：
かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第一式、即ち磁界を表す差分式は、ｘ成分について次の通りとなる。

また、ｙ成分について次の通りとなる。

かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第二式、即ち電界を表す差分式は、次の通りとなる。

以上の各式において、係数αが乗ぜられた項が４次精度３階微分の付加項である。
本ケースの場合、係数αおよびクーラン数を探索するための分散関係式は次式で与えられる。

空間格子がｘ、ｙ方向に等しい間隔（Δｘ＝Δｙ）の場合、この分散関係式に基づいて得られるクーラン数および係数の組み合わせは、次の表の通りである。Δｘ、Δｙが等しくない場合は、ｘ方向、ｙ方向のそれぞれに対して、下表を求めることになる。

既に説明したケースと同様、本ケースにおいても、計算精度を確保または向上しつつ、クーラン数を大きく設定することにより演算の所要時間を短くすることができる。

Ｃ７．３次元４次精度＋２次精度３階微分項：
かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第一式、即ち磁界を表す差分式は、次の通りとなる。

また、Ｍａｘｗｅｌｌ方程式の第二式、即ち電界を表す差分式は、次の通りとなる。

以上の２式において、係数αが乗ぜられた項が４次精度３階微分の付加項である。
本ケースの場合、係数αおよびクーラン数を探索するための分散関係式は次式で与えられる。

空間格子がｘ、ｙ、ｚ方向に等しい間隔（Δｘ＝Δｙ＝Δｚ）の場合、この分散関係式に基づいて得られるクーラン数および係数αの組み合わせは、次の表の通りである。Δｘ、Δｙ、Δｚが等しくない場合は、ｘ方向、ｙ方向、ｚ方向のそれぞれに対して、下表を求めることになる。

既に説明したケースと同様、本ケースにおいても、計算精度を確保または向上しつつ、クーラン数を大きく設定することにより演算の所要時間を短くすることができる。

Ｃ８．２次元４次精度＋２次精度３階微分項：
かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第一式、即ち磁界を表す差分式は、ｘ成分について次の通りとなる。

また、ｙ成分について次の通りとなる。

かかるケースの場合、Ｍａｘｗｅｌｌ方程式の第二式、即ち電界を表す差分式は、次の通りとなる。

以上の各式において、係数αが乗ぜられた項が４次精度３階微分の付加項である。
本ケースの場合、係数αおよびクーラン数を探索するための分散関係式は次式で与えられる。

空間格子がｘ、ｙ方向に等しい間隔（Δｘ＝Δｙ）の場合、この分散関係式に基づいて得られるクーラン数および係数の組み合わせは、次の表の通りである。Δｘ、Δｙが等しくない場合は、ｘ方向、ｙ方向のそれぞれに対して、下表を求めることになる。

既に説明したケースと同様、本ケースにおいても、計算精度を確保または向上しつつ、クーラン数を大きく設定することにより演算の所要時間を短くすることができる。

以上で説明した実施例の解析方法によれば、電磁波または音波の解析を精度よく高速に計算することができる。

音波の伝搬の解析への適用について説明する。
音波の伝搬を記述する簡略化した流体方程式は次式で与えられる。

上式は、Ｍａｘｗｅｌｌ方程式と同じ形式である。従って、実施例で説明した種々の差分式の変数を置き換えることにより、音波の伝搬の解析に適用することができる。
なお、本発明は、本実施例に示した態様に関わらず、その趣旨を変更しない範囲で種々の態様により実現可能である。

本発明は、電磁波または音波の伝搬の解析に利用することができる。

１１ ：入力部
１２ ：算出結果保持部
１３ ：出力部
１４ ：磁界演算部
１５ ：電界演算部
１６ ：係数等設定部
１７ ：設定テーブル

Claims (11)

1. コンピュータによって、ＦＤＴＤ法により電波または音波の伝搬を表す２つの方程式に含まれる２つの状態量を数値計算するＦＤＴＤ解析方法であって、
   前記コンピュータが実行する工程として、
   （ａ）空間内に設定された各格子点における２つの前記状態量その他の初期条件を入力し、前記コンピュータのメモリの所定領域に格納するステップと、
   （ｂ）前記方程式の第一式を表す差分式によって、所定の時間刻みΔｔ後の前記格子点における第１の状態量を計算し、結果を前記コンピュータのメモリの所定領域に格納するステップと、
   （ｃ）前記方程式の第二式を表す差分式によって、所定の時間刻みΔｔ後の前記格子点における第２の状態量を計算し、結果を前記コンピュータのメモリの所定領域に格納するステップと、
   （ｄ）所定の終了条件に至るまで前記ステップ（ｂ）、（ｃ）を繰り返し実行するステップとを備え、
   前記差分式は、ＦＤＴＤ法の計算式に、奇数階微分に係数を乗じた付加項を１つ以上加算した式であるＦＤＴＤ解析方法。
2. 請求項１記載のＦＤＴＤ解析方法であって、
   前記差分式は、ＦＤＴＤ法におけるｎ（ｎは偶数）次精度の計算式に、（ｎ－２ｍ）（ｍはｎ／２より小さい自然数）次以下の精度の（２ｍ＋１）階微分に係数を乗じた付加項を１つ以上加算した式であるＦＤＴＤ解析方法。
3. 請求項１記載のＦＤＴＤ解析方法であって、
   前記付加項に乗ぜられるそれぞれの係数は、クーラン数Ｃと相関を有するＦＤＴＤ解析方法。
4. 請求項３記載のＦＤＴＤ解析方法であって、
   前記ステップ（ａ）において、前記係数またはクーラン数Ｃのいずれかを入力し、
   前記ステップ（ｂ）および前記ステップ（ｃ）の実行前に、前記係数とクーラン数との間の予め設定された相関を利用して、残余の前記係数およびクーラン数を設定するステップを備えるＦＤＴＤ解析方法。
5. 請求項１記載のＦＤＴＤ解析方法であって、
   前記付加項に乗ぜられるそれぞれの係数は、前記伝搬における位相速度の誤差を最小とするという評価基準に基づいて設定されているＦＤＴＤ解析方法。
6. 請求項１記載のＦＤＴＤ解析方法であって、
   前記差分式は、ＦＤＴＤ法における６次精度の計算式に、４次精度の３階微分に係数を乗じた付加項を加算した式であるＦＤＴＤ解析方法。
7. 請求項６記載のＦＤＴＤ解析方法であって、
   前記差分式は、さらに、２次精度の５階微分に係数を乗じた付加項を加算した式であるＦＤＴＤ解析方法。
8. 請求項１記載のＦＤＴＤ解析方法であって、
   前記差分式は、ＦＤＴＤ法における４次精度の計算式に、４次精度の３階微分に係数を乗じた付加項を加算した式であるＦＤＴＤ解析方法。
9. 請求項１記載のＦＤＴＤ解析方法であって、
   前記差分式は、ＦＤＴＤ法における４次精度の計算式に、２次精度の３階微分に係数を乗じた付加項を加算した式であるＦＤＴＤ解析方法。
10. コンピュータによって、ＦＤＴＤ法により電波または音波の伝搬を表す２つの方程式に含まれる２つの状態量を数値計算するためのコンピュータプログラムであって、
    （ａ）空間内に設定された各格子点における２つの前記状態量その他の初期条件を入力し、前記コンピュータのメモリの所定領域に格納する機能と、
    （ｂ）前記方程式の第一式を表す差分式によって、所定の時間刻みΔｔ後の前記格子点における第１の状態量を計算し、結果を前記コンピュータのメモリの所定領域に格納する機能と、
    （ｃ）前記方程式の第二式を表す差分式によって、所定の時間刻みΔｔ後の前記格子点における第２の状態量を計算し、結果を前記コンピュータのメモリの所定領域に格納する機能と、
    （ｄ）所定の終了条件に至るまで前記機能（ｂ）、（ｃ）を繰り返し実行する機能とをコンピュータに実現させることができ、
    前記差分式は、ＦＤＴＤ法における計算式に、奇数階微分に係数を乗じた付加項を１つ以上加算した式であるコンピュータプログラム。
11. ＦＤＴＤ法により電波または音波の伝搬を表す２つの方程式に含まれる２つの状態量を数値計算するためのＦＤＴＤ解析装置であって、
    空間内に設定された各格子点における２つの前記状態量を記憶する算出結果保持部と、
    前記各格子点における前記状態量その他の初期条件を入力し、前記算出結果保持部の所定領域に格納する入力部と、
    前記方程式の第一式を表す差分式によって、所定の時間刻みΔｔ後の前記格子点における第１の状態量を計算し、結果を前記算出結果保持部の所定領域に格納する第１の演算部と、
    前記方程式の第二式を表す差分式によって、所定の時間刻みΔｔ後の前記格子点における第２の状態量を計算し、結果を前記算出結果保持部の所定領域に格納する第２の演算部と、
    所定の終了条件に至るまで前記第１の演算部および第２の演算部による計算を繰り返し実行する制御部とを備え、
    前記差分式は、ＦＤＴＤ法における計算式に、奇数階微分に係数を乗じた付加項を１つ以上加算した式であるＦＤＴＤ解析装置。
    
    
    
    

Images