(57)【要約】【課題】空間微分及び時間発展に関する計算精度を向上させ、電磁界の空間分布の時間変化をより高精度に解析することができる電磁界解析装置を提供する。【解決手段】電界及び磁界並びにその空間微分の初期値を、座標空間に配置された格子点ごとに記憶する初期値記憶部２と、電界及び磁界を格子点間において補間する５次の代数補間式の係数を、初期値に基づいて生成する補間式係数生成部３と、補間式の空間微分からマックスウェルの方程式に基づいて、電界及び磁界の時間微分を生成する時間微分生成部４と、時間微分に基づいて電界及び磁界の時間発展を算出する時間発展算出部５により構成される。【選択図】        図１
2004239784.tif 000002
【特許請求の範囲】
【請求項１】
電界または磁界並びにその空間微分の初期値を、座標空間に配置された格子点ごとに記憶する初期値記憶手段と、
電界または磁界を格子点間において補間する３次以上の補間式の係数を、初期値に基づいて生成する補間式係数生成手段と、
補間式の空間微分からマックスウェルの方程式に基づいて、電界または磁界の時間微分を生成する時間微分生成手段と、
時間微分に基づいて電界または磁界の時間発展を算出する時間発展算出手段を備えたことを特徴とする電磁界解析装置。
【請求項２】
上記補間式係数生成手段は、補間式として５次の代数式の係数を生成し、代数式が座標軸方向に隣接する３つの格子点に対応していることを特徴とする請求項１に記載の電磁界解析装置。
【請求項３】
上記初期値記憶手段は、電界及び磁界と、その空間微分の各初期値を記憶し、
上記補間式係数生成手段は、電界及び磁界の各補間式の係数を生成することを特徴とする請求項１または２に記載の電磁界解析装置。
【請求項４】
上記初期値記憶手段は、電界及びその時間微分と、電界の空間微分及びその時間微分の各初期値を記憶し、
上記補間式係数生成手段は、電界及びその時間微分のそれぞれに対応する補間式の係数を生成することを特徴とする請求項１または２に記載の電磁界解析装置。
【請求項５】
電界または磁界並びにその空間微分の初期値を、座標空間に配置された格子点ごとに記憶する初期値記憶ステップと、
電界または磁界を格子点間において補間する３次以上の補間式の係数を、初期値に基づいて生成する補間式係数生成ステップと、
補間式の空間微分からマックスウェルの方程式に基づいて、電界または磁界の時間微分を生成する時間微分生成ステップと、
時間微分に基づいて電界または磁界の時間発展を算出する時間発展算出ステップからなることを特徴とする電磁界解析方法。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
本発明は、電磁界解析装置及び電磁界解析方法に係り、より詳しくは、電磁界の空間分布の時間変化を解析する電磁界解析装置の改良に関する。
【０００２】
【従来の技術】
従来の電磁界解析方法としては、例えば、特許文献１に記載のものがある。この特許文献１に記載の電磁界解析方法は、有限差分時間領域（Ｆｉｎｉｔｅ  Ｄｉｆｆｅｒｅｎｃｅ  Ｔｉｍｅ  Ｄｏｍａｉｎ：ＦＤＴＤ）法を用いて解析対象物が発生する時間領域の近傍電磁界を求め、近傍電磁界に等価定理を適用して解析対象物からの遠方放射界を計算することにおいて、近傍電磁界をフーリエ変換により周波数領域に変換する際に、サンプリング時間を所望の最大周波数成分までが得られるように間引くことにより、遠方放射界の計算に必要なメモリ量が低減するように構成される。
【０００３】
また、流体解析におけるＩＤＯ（Ｉｎｔｅｒｐｏｌａｔｅｄ  Ｄｉｆｆｅｒｅｎｔｉａｌ  Ｏｐｅｒａｔｏｒ）法については、非特許文献１に記載されている。
【０００４】
【特許文献１】
特開平１１−２９５３６５号公報
【非特許文献１】
「インターポレイテッド・ディファレンシャル・オペレータ・スキーム・フォー・ソルビング・パーシャル・ディファレンシャル・イクエーションズ（Ｉｎｔｅｒｐｏｌａｔｅｄ  Ｄｉｆｆｅｒｅｎｔｉａｌ  Ｏｐｅｒａｔｏｒ  ｓｃｈｅｍｅ  ｆｏｒ  ｓｏｌｖｉｎｇ  ｐａｒｔｉａｌ  ｄｉｆｆｅｒｅｎｔｉａｌ  ｅｑｕａｔｉｏｎｓ）」、コンピュータ・フィジックス・コミュニケーションズ（Ｃｏｍｐｕｔｅｒ  Ｐｈｙｓｉｃｓ  Ｃｏｍｍｕｎｉｃａｔｉｏｎｓ）、１９９７年、１０２、ｐ．１３２−１４６
【０００５】
【発明が解決しようとする課題】
上記したような従来のＦＤＴＤ法による電磁界の解析は、電界Ｅ及び磁界Ｈという２種類の物理量に対して立てられたマックスウェルの方程式に基づいて、時間微分を空間微分に置き換えて時間発展を求めようとするものであり、空間微分に差分法を用いている。このため、座標空間に配置された格子点間を直線により補間することに相当し、空間微分及び時間発展に関して１次の近似精度までしか解析することができない。
【０００６】
一方、従来のＩＤＯ法による解析は、１つの物理量に対して立てられた流体方程式または波動方程式を解こうとするものである。例えば、ｆを物理量、ｔを時間、ｘを空間座標、ｃを定数、ｇを補助的な関数として、∂ｆ／∂ｔ＋ｃ∂ｆ／∂ｘ＝ｇという形で代表される方程式をｆについて解析している。従って、２種類の物理量に対するマックスウェルの方程式を解析するものではない。
【０００７】
また、マックスウェルの方程式は、電界Ｅの空間微分及び磁界Ｈの時間微分が結合した形の式と、Ｈの空間微分及びＥの時間微分が結合した形の式との連立偏微分方程式である。流体方程式も連立偏微分方程式であるが、上記した式において、ｆとして密度ρをとり、ｇの中に速度ｕが入った式と、ｆとして速度ｕをとり、ｇの中に密度ρが入った式との連立という形態であり、マックスウェルの方程式とは本質的に異なっている。
【０００８】
本発明は、上記事情に鑑みてなされたものであり、空間微分及び時間発展に関する計算精度を向上させ、電磁界の空間分布の時間変化をより高精度に解析することができる電磁界解析装置及び電磁界解析方法を提供することを目的としている。
【０００９】
【課題を解決するための手段】
本発明による電磁界解析装置は、電界または磁界並びにその空間微分の初期値を、座標空間に配置された格子点ごとに記憶する初期値記憶手段と、電界または磁界を格子点間において補間する３次以上の補間式の係数を、初期値に基づいて生成する補間式係数生成手段と、補間式の空間微分からマックスウェルの方程式に基づいて、電界または磁界の時間微分を生成する時間微分生成手段と、時間微分に基づいて電界または磁界の時間発展を算出する時間発展算出手段により構成される。
【００１０】
この様な構成によれば、３次以上の補間式を用いるので、空間微分及び時間発展に関する計算精度が向上し、電磁界の空間分布の時間変化をより高精度に解析することができる。
【００１１】
本発明による電磁界解析方法は、電界または磁界並びにその空間微分の初期値を、座標空間に配置された格子点ごとに記憶する初期値記憶ステップと、電界または磁界を格子点間において補間する３次以上の補間式の係数を、初期値に基づいて生成する補間式係数生成ステップと、補間式の空間微分からマックスウェルの方程式に基づいて、電界または磁界の時間微分を生成する時間微分生成ステップと、時間微分に基づいて電界または磁界の時間発展を算出する時間発展算出ステップからなる。
【００１２】
【発明の実施の形態】
実施の形態１．
図１は、本発明の実施の形態１による電磁界解析装置の一構成例を示したブロック図である。本実施の形態の電磁界解析装置１は、電磁回路の解析などに用いられ、電磁界の解析手法は、電磁界の空間分布及びその時間変化を決める基本方程式であるマックスウェルの方程式を解くために、ＩＤＯ法及びＦＤＴＤ法の組み合わせを基礎としている。そこで、本実施の形態における電磁界の解析手法を理解するために、まず、ＩＤＯ法及びＦＤＴＤ法について説明する。
【００１３】
図２は、座標空間に配置された格子点ごとの物理量に対して物理量を格子点間において補間する代数曲線を示した説明図である。流体解析におけるＩＤＯ法は、物理現象を記述する微分方程式を解くために考案されたものである。そのＩＤＯ法では、まず座標空間に格子を張り、格子点上に物理量と、その物理量の空間微分量が初期値として配置される。そして、格子点間の物理量を３次または５次の代数曲線により補間する。
【００１４】
ここでは、５次の代数曲線により補間する場合について説明し、空間の次元を１次元とする。座標軸（ｘ軸）上のＮ＋１個の格子点のうち、ｉ番目の格子点に着目し、この格子点を挟んでｘ軸方向に隣接する３つの格子点間を補間する代数曲線は、次式（１）の補間式Ｆｉ（ｘ）により表される。
【数１】
2004239784.tif 000003
【００１５】
式（１）において、ｉ＝１〜Ｎ−１としてｉ番目の格子点における物理量ｆの初期値をｆｉとし、その空間微分（ｘ方向微分）ｆｘの初期値をｆｘｉとしている。補間式Ｆｉ（ｘ）の係数ａｉ〜ｄｉは、格子点の間隔をΔｘとすると、次式（２ａ）〜（２ｄ）により表され、補間式Ｆｉ（ｘ）がこの区間内において微分方程式を満たす解となるための条件Ｆｉ（Δｘ）＝ｆｉ＋１，Ｆｉ（−Δｘ）＝ｆｉ−１，Ｆｘｉ（Δｘ）＝ｆｘｉ＋１，Ｆｘｉ（−Δｘ）＝ｆｘｉ−１から求められる。
【数２】
2004239784.tif 000004
【００１６】
また、式（１）から空間の１階微分から５階微分は、それぞれ次式（３）〜（７）により表され、結局、物理量ｆｉ及びその空間微分量ｆｘｉが微分方程式を用いて格子点ごとに求められていれば、補間式Ｆｉ（ｘ）は微分方程式の解となっている。
【数３】
2004239784.tif 000005
【００１７】
次に初期時刻ｔからΔｔ後の時間発展を求めるために、物理量ｆ及びその空間微分量ｆｘをそれぞれテーラー展開する。ここで、物理量ｆはΔｔについて５次の項までとり、空間微分量ｆｘは４次の項までとれば十分な近似精度が得られる。この様にして求められる時間ステップΔｔ後のｆ及びｆｘは、それぞれ次式（８）及び（９）により表される。
【数４】
2004239784.tif 000006
【００１８】
これらの式（８）及び（９）における時間の１階微分から５階微分までの各項が求められると、物理量ｆ及び空間微分ｆｘの時間発展をそれぞれ求めることができる。ＩＤＯ法では、微分方程式に基づいて時間微分を空間微分に変換することにより、これらの時間微分の各項を求めている。
【００１９】
微分方程式がｋを定数として、次式（１０）により表される場合、ｆの時間の１階微分は、この式（１０）を用いて空間微分に変換される。
【数５】
2004239784.tif 000007
【００２０】
補間式Ｆｉ（ｘ）は微分方程式（１０）を満たすので、変換された空間微分は式（３）を用いてその値を求めることができる。時間の２階微分については、式（１０）の両辺を時間ｔで微分して、次式（１１）を求める。
【数６】
2004239784.tif 000008
【００２１】
そして、再度式（１０）を用いて時間微分を空間微分に変換すれば、時間の２階微分は、空間微分だけで表され、式（４）を用いてその値を求めることができる。時間の３〜５階微分についても同様にして、式（５）〜（７）を用いて値を求めることができる。
【００２２】
また、ｆｘの時間の１階微分については、式（１０）の両辺を空間ｘで微分して次式（１２）を求めれば、この式（１２）を用いて空間微分に変換される。従って、式（４）を用いてその値を求めることができる。
【数７】
2004239784.tif 000009
【００２３】
ｆｘの時間の２階微分については、式（１２）の両辺を時間ｔで微分して、次式（１３）を求める。
【数８】
2004239784.tif 000010
【００２４】
そして、再度式（１２）を用いて時間微分を空間微分に変換すれば、ｆｘの時間の２階微分は、空間微分だけで表され、式（５）を用いて値を求めることができる。ｆｘの時間の３及び４階微分も同様にして、式（６）及び（７）を用いて求められる。
【００２５】
この様にして、時間微分項を求めることにより、式（８）及び（９）を用いて物理量ｆｉ及び空間微分ｆｘｉの時間発展を求めることができる。時間発展が求められれば、これらの時間発展を初期値として、再び式（２ａ）〜（２ｄ）に戻って係数ａｉ〜ｄｉを求めるというようにして、次々と微分方程式（１０）を満たす物理量ｆの時間発展を求めていくことができる。
【００２６】
図３は、座標空間における格子点ごとの電界及び磁界の配置を示した説明図である。ＦＤＴＤ法は、差分法を用いてマックスウェルの方程式を直接的に解く方法である。そのＦＤＴＤ法では、ＩＤＯ法と同様にまず座標空間に格子を張って、そこに電界及び磁界の座標成分が初期値として配置される。
【００２７】
座標成分の配置については、空間差分の近似精度を上げるため、同一の格子点には電界及び磁界を配置せずに、電界が配置された格子点間の中央の位置に磁界が配置され、磁界が配置された格子点間の中央の位置に電界が配置される。つまり、電界及び磁界が座標軸方向に交互に配置される。また、時間ステップ幅をΔｔとすると、電界及び磁界を求める時刻は、Δｔ／２だけずらして求められる。
【００２８】
ｘ方向及びｙ方向には空間的に変化しない平面波の場合、電界のｘ方向成分をＥｘとし、磁界のｙ方向成分をＨｙとすると、初期時刻ｔに対してΔｔ後の時間発展は、マックスウェルの方程式に基づいて、それぞれ次式（１４）及び（１５）により求められる。
【数９】
2004239784.tif 000011
【００２９】
式（１４）及び（１５）において、εは誘電率であり、μは透磁率である。計算の手順としては、式（１４）により電界Ｅｘｉの時間発展を求め、式（１５）により磁界Ｈｙｉの時間発展を求める。そして、これらを初期値として再び式（１４）及び（１５）により電界Ｅｘｉ及び磁界Ｈｙｉの時間発展を求めていくというようにして行われる。
【００３０】
この様に、ＦＤＴＤ法では、空間微分に差分を用いているので、電界及び磁界を格子点間において直線により補間することに相当する。このため、時間に関して、テーラー展開の１次の時間微分項までしか得られず、ＩＤＯ法に比べて精度が低い。
【００３１】
次に、本実施の形態の電磁界解析装置１における解析手法について説明する。この電磁界解析装置１は、初期値記憶部２、補間式係数生成部３、時間微分生成部４、時間発展算出部５、ディスプレイ６及び補助記憶部７により構成される。
【００３２】
初期値記憶部２は、電界及び磁界の初期値と、電界及び磁界の空間微分の初期値を格納する記憶手段であり、これらの初期値は、座標空間に配置された格子点ごとに予め記憶されている。
【００３３】
補間式係数生成部３は、電界及び磁界を格子点間において補間する補間式の係数を、初期値に基づいてそれぞれ算出している。各補間式は、３次以上の代数式、例えば、５次の代数式により表され、座標軸方向に隣接する３つの格子点に対応している。生成した補間式の係数は、格子点ごとに時間微分生成部４へ出力される。
【００３４】
時間微分生成部４は、電界及び磁界の時間微分と、電界及び磁界の空間微分に対する時間微分をそれぞれ算出している。これらの時間微分は、補間式の空間微分からマックスウェルの方程式に基づいて算出され、電界及び磁界については時間に関して５階までの微分量が生成され、電界及び磁界の空間微分については時間に関して４階までの微分量が生成される。生成した時間微分は、格子点ごとに時間発展算出部５へ出力される。
【００３５】
時間発展算出部５は、時間微分量から、電界及び磁界と、それらの空間微分の時間発展をそれぞれ算出している。これらの時間発展は、初期時刻から微小時間後についてのテーラー展開に基づいて格子点ごとに算出される。電界及び磁界の時間発展については、微小時間について５次の項までテーラー展開され、電界及び磁界の空間微分の時間発展については、４次の項までテーラー展開される。それぞれ必要に応じて、上記したものより低次の項までのテーラー展開であっても良い。
【００３６】
算出した時間発展は、算出結果としてディスプレイ６へ出力されて画面表示されるとともに、補助記憶部７に記憶される。時間発展算出部５により格子点ごとの電界及び磁界の時間発展と、電界及び磁界の空間微分についての時間発展が算出されると、これらの時間発展は補間式係数生成部３へ出力され、これらを新たな初期値として補間式の係数が生成される。
【００３７】
図４は、座標空間における格子点ごとの電界及び磁界並びにそれらの空間微分の配置の一例を示した説明図である。マックスウェルの方程式は、電界ベクトルをＥとし、磁界ベクトルをＨとして、次式（１６）〜（１９）により表される。
【数１０】
2004239784.tif 000012
【００３８】
式（１７）及び（１８）において、σは導電率であり、ρは電荷密度である。有限要素法のような計算手法では、通常、電磁波の伝搬を解析する場合に、上式（１６）〜（１９）を演算し、電界ベクトルＥまたは磁界ベクトルＨのいずれか一方のみについての波動方程式もしくはヘルムホルツ方程式の形式に変形して計算が行われるが、本実施の形態では、マックスウェルの方程式（１６）〜（１９）がそのままの形で計算に用いられる。
【００３９】
以下、空間座標にｘ，ｙ，ｚ直交座標を用い、ｘ方向及びｙ方向には電磁界成分が空間的に変化しない平面波の伝搬現象を解析する場合について説明する。そして、空間中には電荷が存在しないものとして、ρ＝０とする。また、導体外部における伝搬現象を解析するものとして、σ＝０とする。さらに、誘電率ε及び透磁率μはそれぞれ一定で空間的及び時間的に変化しないものとする。これらの条件は説明を簡単にするためのものであり、本発明を制限するものではない。
【００４０】
このとき、電界ベクトルＥのｘ方向成分をＥｘ、ｙ方向成分をＥｙ、ｚ方向成分をＥｚとし、磁界ベクトルＨのｘ方向成分をＨｘ、ｙ方向成分をＨｙ、ｚ方向成分をＨｚとすると、式（１６）及び（１７）は、次式（２０ａ）〜（２１ｃ）によりそれぞれ表される。
【数１１】
2004239784.tif 000013
【００４１】
式（１８）及び（１９）については、電界ベクトルＥ及び磁界ベクトルＨの初期値における空間分布がこれらの式を満たせば、その後の時間発展の演算には必要ない。式（２０ｃ）からＥｚ＝０であり、式（２１ｃ）からＨｚ＝０であることがわかる。また、式（２０ａ）〜（２１ｃ）は、式（２０ａ）及び（２１ｂ）の組と、式（２０ｂ）及び（２１ａ）の組とに分離できることがわかる。
【００４２】
そこで、座標空間のｚ軸方向にＮ＋１個の格子点を配置し、これらの格子点上に電界ベクトルＥ及び磁界ベクトルＨの各座標成分Ｅｘ，Ｅｙ，Ｈｙ，Ｈｘの初期値と、電界及び磁界の空間微分（ｚ方向微分）の初期値をそれぞれ配置する。このとき、電界及び磁界の座標成分と、それらの空間微分とはすべて同一の格子点上に配置される。
【００４３】
例えば、ｉ番目の格子点に着目すると、この格子点には、初期値として、Ｅｘｉ及びその空間微分量Ｅｘｚｉと、Ｅｙｉ及びその空間微分量Ｅｙｚｉと、Ｈｙｉ及びその空間微分量Ｈｙｚｉと、Ｈｘｉ及びその空間微分量Ｈｘｚｉが配置されている（ｉ＝０〜Ｎ）。
【００４４】
格子点の間隔は、必ずしも等間隔である必要はないが、ここでは、一定のΔｚという間隔をとるものとする。そして、以下では、式（２０ａ）及び（２１ｂ）の組により表される電磁界現象の解析を例として説明する。式（２０ｂ）及び（２１ａ）の組についても同様である。
【００４５】
図５及び図６のステップＳ１〜Ｓ１３は、図１の電磁界解析装置における電磁界解析の動作の一例を示したフローチャートである。まず、電磁界の解析処理の準備として、格子点の間隔Δｚと、時間ステップ幅Δｔと、解析処理の終了時刻ｔＥＮＤと、ｉ＝０〜Ｎとして格子点上の電界Ｅ及び磁界Ｈの座標成分の初期値Ｅｘｉ及びＨｙｉと、電界Ｅ及び磁界Ｈの空間微分の初期値Ｅｘｚｉ及びＨｙｚｉを、初期値記憶部２に予め記憶させておく（ステップＳ１）。
【００４６】
また、補間式係数生成部３内において、Δｚと、初期値Ｅｘｉ，Ｈｙｉ，Ｅｘｚｉ，Ｈｙｚｉ、或いは、これらの更新された値ＥｘｉＮＥＷ，ＨｙｉＮＥＷ，ＥｘｚｉＮＥＷ，ＨｙｚｉＮＥＷに対応するデータを記憶するメモリ領域を確保する。
【００４７】
ここで、時間ステップ幅Δｔは、空間中を伝搬する光の速度をｃとして、次式（２２）により求められる。
【数１２】
2004239784.tif 000014
【００４８】
式（２２）において、光の速度ｃは、誘電率ε及び透磁率μから次式（２３）により求められ、ＣＦＬはクーラン条件と呼ばれ、計算精度を上げるために１以下の定数とされる。
【数１３】
2004239784.tif 000015
【００４９】
以上のような計算の準備と初期化が終われば、初期時刻をｔ＝０として、解析を開始する（ステップＳ２）。まずｉ＝１としてｉ番目の格子点に着目する（ステップＳ３）。このとき、必要であれば、電磁界の励振条件を算出する。
【００５０】
次に、補間式係数生成部３は、ｉ番目の格子点を挟んでｚ軸方向に隣接する３つの格子点間においてＥｘｉ及びＨｙｉをそれぞれ補間する補間式の係数を算出する（ステップＳ４）。補間式は、次式（２４）及び（２５）によりそれぞれ定義される。
【数１４】
2004239784.tif 000016
【００５１】
式（２４）及び（２５）において、係数ＡＥ〜ＤＥ及びＡＨ〜ＤＨは、補間式がこの区間内においてマックスウェルの方程式の解となるための条件から求められる。例えば、係数ＡＥ〜ＤＥは、ｉ＝１〜Ｎ−１として、次式（２６ａ）〜（２６ｄ）により表される連立方程式を解くことによってそれぞれ求められる。
【数１５】
2004239784.tif 000017
【００５２】
係数ＡＥ〜ＤＥは、式（２６ａ）〜（２６ｄ）に基づいて、次式（２７ａ）〜（２７ｄ）により表される。
【数１６】
2004239784.tif 000018
【００５３】
係数ＡＨ〜ＤＨについても、式（２６ａ）〜（２６ｄ）と同様の条件からそれぞれ求められ、次式（２８ａ）〜（２８ｄ）により表される。
【数１７】
2004239784.tif 000019
【００５４】
補間式の係数ＡＥ〜ＤＥ及びＡＨ〜ＤＨは、これらの式（２７ａ）〜（２８ｄ）に基づいて算出される。
【００５５】
時間微分生成部４は、算出された補間式の係数ＡＥ〜ＤＥ及びＡＨ〜ＤＨに基づいて、Ｅｘ，Ｅｘｚ，Ｈｙ及びＨｙｚの時間微分量をそれぞれ算出する（ステップＳ５）。これらの時間微分量は、以下のようにして求められる。
【００５６】
Ｅｘの時間に関する１階微分については、まず、式（２０ａ）を次式（２９）のように変形する。
【数１８】
2004239784.tif 000020
【００５７】
式（２５）を空間微分（ｚで微分）することにより、式（２９）から次式（３０）が求められる。
【数１９】
2004239784.tif 000021
【００５８】
また、Ｅｘの時間の２階微分について、式（２９）を時間微分して式（２１ｂ）を用いることにより、次式（３１）が求められる。
【数２０】
2004239784.tif 000022
【００５９】
従って、式（２４）を２階空間微分することにより、式（３１）から次式（３２）が得られる。
【数２１】
2004239784.tif 000023
【００６０】
さらに、Ｅｘの時間の３階微分について、式（３１）を時間微分して式（２９）を用いることにより、次式（３３）が求められる。
【数２２】
2004239784.tif 000024
【００６１】
従って、式（２５）を３階空間微分することにより、式（３３）から次式（３４）が得られる。
【数２３】
2004239784.tif 000025
【００６２】
Ｅｘの時間の４階微分についても、式（３３）を時間微分して式（２１ｂ）を用いれば、次式（３５）が求められる。
【数２４】
2004239784.tif 000026
【００６３】
従って、式（２４）を４階空間微分することにより、式（３５）から次式（３６）が得られる。
【数２５】
2004239784.tif 000027
【００６４】
Ｅｘの時間の５階微分についても、式（３５）を時間微分して式（２９）を用いれば、次式（３７）が求められる。
【数２６】
2004239784.tif 000028
【００６５】
従って、式（２５）を５階空間微分すれば、式（３７）から次式（３８）が得られる。
【数２７】
2004239784.tif 000029
【００６６】
Ｅｘｚの時間に関する１〜４階微分については、式（３０），（３２），（３４），（３６）を空間微分することにより、それぞれ次式（３９）〜（４２）によって表される。
【数２８】
2004239784.tif 000030
【００６７】
結局、ｉ番目の格子点におけるＥｘ及びＥｘｚの時間微分量は、ｚ＝０として式（３０），（３２），（３４），（３６），（３８），（３９）〜（４２）から得られる次式（４３）及び（４４）によりそれぞれ算出される。
【数２９】
2004239784.tif 000031
【００６８】
式（４３）において、Ｅｘｔ〜Ｅｘｔｔｔｔｔは、それぞれＥｘの時間に関する１〜５階微分量を表し、式（４４）において、Ｅｘｚｔ〜Ｅｘｚｔｔｔｔは、それぞれＥｘｚの時間に関する１〜４階微分量を表す。
【００６９】
Ｈｙ及びＨｙｚの時間微分量についても同様にして、次式（４５）及び（４６）により算出される。
【数３０】
2004239784.tif 000032
【００７０】
式（４５）において、Ｈｙｔ〜Ｈｙｔｔｔｔｔは、それぞれＨｙの時間に関する１〜５階微分量を表し、式（４６）において、Ｈｙｚｔ〜Ｈｙｚｔｔｔｔは、それぞれＨｙｚの時間に関する１〜４階微分量を表す。
【００７１】
次に、時間発展算出部５は、算出された時間微分量に基づいて、Ｅｘ，Ｅｘｚ，Ｈｙ及びＨｙｚの時間発展をそれぞれ算出する（ステップＳ６）。初期時刻ｔからΔｔ後の時間発展は、以下のようにして求められる。
【００７２】
電界のｘ方向成分Ｅｘ及び磁界のｙ方向成分Ｈｙについては、時間に関して５次の近似精度までテーラー展開し、空間微分Ｅｘｚ及びＨｙｚについては、４次の精度までテーラー展開することにより、次式（４７）〜（５０）が得られる。
【数３１】
2004239784.tif 000033
【００７３】
式（４７）〜（５０）における時間微分量に式（４３）〜（４６）をそれぞれ用いることにより、ｉ番目の格子点における時刻ｔ＋Δｔでの値が、次式（５１）〜（５４）から求められる。
【数３２】
2004239784.tif 000034
【００７４】
式（５１）〜（５４）によって時間発展が算出され、算出結果は、時間発展算出部５により補助記憶部７に記憶される。
【００７５】
以上のようにして、ｉ番目の格子点における時間発展が求められると、ｉ＝ｉ＋１とし（ステップＳ７）、次の格子点における時間発展を求めるために、ステップＳ４〜７が繰り返される。そして、ｉ＝Ｎ−１となるまで繰り返した後は、境界条件が時間発展算出部５により処理される（ステップＳ８，Ｓ９）。
【００７６】
境界条件の処理は、ｉ＝０及びｉ＝ＮにおけるＥｘｉ，Ｅｘｚｉ，Ｈｙｉ及びＨｙｚｉの時間発展を定めることにより行われ、例えば、それぞれ０であると定められる。
【００７７】
次に、時間発展算出部５は、生成した全ての格子点（ｉ＝０〜Ｎ）における時間発展を、Ｅｘｉ，Ｅｘｚｉ，Ｈｙｉ及びＨｙｚｉの更新値として、補助記憶部７に記憶するとともに、ディスプレイ６に出力して画面表示させる（ステップＳ１０）。
【００７８】
また、これらの更新値は、時間発展算出部５により補間式係数生成部３に出力され（ステップＳ１１）、補間式係数生成部３は、ｔ＝ｔ＋Δｔとして（ステップＳ１２）、更新値に基づいて補間式の係数を算出する。そして、ステップＳ３〜Ｓ１２が繰り返される。
【００７９】
この様にして、次々に時間発展が求められ、時刻ｔがｔＥＮＤ以上になると（ステップＳ１３）、この解析処理は終了される。
【００８０】
図７は、図１の電磁界解析装置による平面矩形波の伝搬現象の解析結果の一例を示した状態図であり、図７（ａ）は、初期時刻（ｔ＝０）における電界Ｅｘを表し、図７（ｂ）は、３０ステップ後（ｔ＝３０Δｔ）、図７（ｃ）は、１００ステップ後（ｔ＝１００Δｔ）、図７（ｄ）は、２００ステップ後（ｔ＝２００Δｔ）の電界Ｅｘをそれぞれ表している。
【００８１】
電界Ｅｘの初期値を、ｚ軸上に配置された格子点うち、中央部の４０〜６０番目の格子点において１とし、その他の格子点において０としている。この様な平面矩形波は、解析処理の時間ステップに関して３０ステップ後には、初期値の１／２の強度の矩形波がもとの矩形波から左右に分離を開始し、１００ステップ後には、ほぼ分離を終えて左右方向にそれぞれ伝搬している。そして、２００ステップ後には、完全に分離し、ほぼ矩形を保った状態のまま伝搬していることがわかる。
【００８２】
図８は、従来のＦＤＴＤ法による平面矩形波の伝搬現象の解析結果の一例を示した状態図であり、図８（ａ）は、初期時刻（ｔ＝０）における電界Ｅｘを表し、図８（ｂ）は、３０ステップ後、図８（ｃ）は、１００ステップ後の電界Ｅｘをそれぞれ表している。初期条件は、図７におけるのと同様である。３０ステップ後において既に矩形波の両端部が発散しており、１００ステップ後には、解析不能となっている。従って、ＦＤＴＤ法では、矩形波の伝搬の様子が正しく解析できないことがわかる。
【００８３】
本実施の形態によれば、空間微分及び時間発展に関する計算精度が向上し、電磁界の空間分布の時間変化を高精度に解析することができる。具体的には、５次の補間式を格子点間における補間に用いたことにより、電界及び磁界の空間微分及び時間発展に関して５次の近似精度を確保することができる。
【００８４】
また、従来のＦＤＴＤ法では、電界及び磁界を格子点ごとに交互に配置しなければならないので、例えば、電界が配置された格子点（位置）における磁界を求めるためには、補間計算をしなければならなかった。本実施の形態では、同一の格子点に電界及び磁界を配置するので、その様な補間計算は不要である。
【００８５】
また、従来のＦＤＴＤ法では、時間に関して中心差分になるように電界及び磁界を互いにΔｔ／２だけずれた時刻で求める。このため、例えば、電界を求めた時刻における磁界を得るためには、時間的な補間計算をしなければならなかった。本実施の形態では、同一の時刻において電界及び磁界を求めるので、その様な補間計算も不要である。
【００８６】
実施の形態２．
図９は、本発明の実施の形態２による電磁界解析装置における格子点ごとの電界及び磁界並びにそれらの空間微分の配置の一例を示した説明図である。実施の形態１では、格子点が座標空間のｚ軸上に等間隔で配置されていたが、本実施の形態では、格子点の間隔が格子点ごとに異なっている。
【００８７】
ｉ番目の格子点に着目して、ｉ−１番目とｉ番目の格子点間の間隔をΔｚｍとし、ｉ番目とｉ＋１番目の格子点間の間隔をΔｚｐとする。このとき、例えば、Ｅｘｉを補間する補間式（２４）の係数ＡＥ〜ＤＥは、次式（５５ａ）〜（５５ｄ）により表される４元連立方程式を解くことによって求められる。
【数３３】
2004239784.tif 000035
【００８８】
式（５５ａ）〜（５５ｄ）は、次式（５６）により表される。
【数３４】
2004239784.tif 000036
【００８９】
この４元連立一次方程式（５６）は、例えば、「フォートラン（Ｆｏｒｔｒａｎ）７７による数値計算ソフトウェア」（丸善株式会社出版）に紹介されているガウス法を用いて解かれ、係数ＡＥ〜ＤＥが補間式係数生成部３によって数値計算される。その他の構成は、実施の形態１と同様である。
【００９０】
本実施の形態によれば、格子点の間隔が等間隔でない場合であっても、電磁界の時間発展を高精度で解析することができる。
【００９１】
実施の形態３．
図１０のステップＳ２１〜Ｓ３２は、本発明の実施の形態３による電磁界解析装置の動作の一例を示したフローチャートである。実施の形態１では、ｘ及びｙ方向には電磁界成分が変化しない平面波の伝搬現象を解析していたが、本実施の形態では、より一般的な場合の解析を行っている。
【００９２】
σ＝０として、マックスウェルの方程式（１６）及び（１７）は、次式（５７ａ）〜（５８ｃ）により表される。
【数３５】
2004239784.tif 000037
【００９３】
座標空間に配置される格子点については、ｘ軸方向にＮｘ＋１個、ｙ軸方向にＮｙ＋１個、ｚ軸方向にＮｚ＋１個の格子点がそれぞれ配置される。格子点の間隔は、必ずしも等間隔である必要はないが、ここでは、それぞれの座標軸に沿って一定であるとし、Δｘ，Δｙ及びΔｚとする。
【００９４】
この様に格子点を配置したとき、第ｉ，ｊ，ｋ番目の格子点に着目して、この格子点上に配置される電界のｘ方向成分をＥｘｉｊｋとし、そのｘ方向微分量をＥｘｘｉｊｋ、ｙ方向微分量をＥｘｙｉｊｋ、ｚ方向微分量をＥｘｚｉｊｋとする。電界のｙ方向及びｚ方向成分についても同様であるとする。また、当該格子点上に配置される磁界のｘ方向成分をＨｘｉｊｋとし、そのｘ方向微分量をＨｘｘｉｊｋ、ｙ方向微分量をＨｘｙｉｊｋ、ｚ方向微分量をＨｘｚｉｊｋとする。磁界のｙ方向及びｚ方向成分についても同様であるとする。
【００９５】
また、電磁界が伝搬する媒質に固有の定数である誘電率ε、透磁率μもしくは導電率σが空間的に変化する場合は、これらを電界及び磁界と同様に各格子点上に配置し、初期値記憶部２に記憶する。
【００９６】
電磁界解析の処理手順において、ステップＳ２１及びＳ２２は、図５のステップＳ１及びＳ２と同様である。ステップＳ２３において、ｉ＝１，ｊ＝１及びｋ＝１とし、ステップＳ２４において、第ｉ，ｊ，ｋ番目の格子点に着目して、補間式係数生成部３により補間式の係数が算出される。
【００９７】
ここで、補間式は、第ｉ，ｊ，ｋ番目の格子点を挟んで各座標軸方向に隣接する３つの格子点に対応しており、例えば、Ｅｘに関しては、次式（５９）により表される。
【数３６】
2004239784.tif 000038
【００９８】
他の電界成分及び磁界成分の補間式も同様に表される。式（５９）のｘ方向微分、ｙ方向微分及びｚ方向微分は、次式（６０ａ）〜（６０ｃ）によりそれぞれ表される。
【数３７】
2004239784.tif 000039
【００９９】
補間式（５９）の係数は、次式（６１ａ）〜（６３ｄ）により表される１２元連立一次方程式を解くことにより求められる。
【数３８】
2004239784.tif 000040
【０１００】
例えば、式（６３ａ）〜（６３ｄ）を書き下すと、次式（６４ａ）〜（６４ｄ）がそれぞれ得られる。
【数３９】
2004239784.tif 000041
【０１０１】
これらの式（６４ａ）〜（６４ｄ）は、式（２６ａ）〜（２６ｄ）を書き下したものと同じ形式であるので、係数Ａｚ〜Ｄｚは、式（２７ａ）〜（２７ｄ）と同様の形式により算出することができる。他の係数Ａｘ〜Ｄｘ及びＡｙ〜Ｄｙについても、式（６１ａ）〜（６２ｄ）から同様にして算出される。電界の他の座標成分についても同様であり、磁界についても同様である。
【０１０２】
ステップＳ２５及びＳ２６は、図５のステップＳ５及びＳ６における処理手順と同様にして行われ、ステップＳ２７において、ｋ＝ｋ＋１とし、ｋがＮｚ以上となるまで、ステップＳ２４〜Ｓ２７が繰り返される。
【０１０３】
ｋがＮｚ以上となると、ｊ＝ｊ＋１，ｋ＝１とし（ステップＳ２８，Ｓ２９）、ｊがＮｙ以上となるまで、ステップＳ２４〜Ｓ２９が繰り返される。そして、ｊがＮｙ以上となると、ｉ＝ｉ＋１，ｊ＝１，ｋ＝１とし（ステップＳ３０，Ｓ３１）、ｉがＮｘ以上となるまで、ステップＳ２４〜Ｓ３１が繰り返される。
【０１０４】
この様にして、境界（ｉ＝０，ｉ＝Ｎｘ，ｊ＝０，ｊ＝Ｎｙ，ｋ＝０、または、ｋ＝Ｎｚ）を除いて座標空間内の全ての格子点について電磁界の時間発展が算出される。その他の処理手順は、実施の形態１と同様にして行われる。
【０１０５】
本実施の形態によれば、３次元のより一般的な電磁界現象に対する精度の高い解析を行うことができる。
【０１０６】
実施の形態４．
図１１のステップＳ４１〜Ｓ５０は、本発明の実施の形態４による電磁界解析装置の動作の一例を示したフローチャートであり、微分演算の処理手順を表している。実施の形態３では、予め求められた関係式に基づいて、時間微分量が算出されるが、本実施の形態では、微分演算が記号的な演算により行われて時間微分量が算出される。
【０１０７】
まず、数式の表現方法について説明する。９つの数値を１組にして「項」を表し、その１組を配列形式により表す。配列形式をｆとすれば、ｎ＝０〜８として、ｎ番目の配列の要素は、ｆ［ｎ］により表される。ｆ［０］には、項の符号「＋」または「−」に対応する数値を割り当て、例えば、「＋」のときは、１が割り当てられ、「−」のときは、−１が割り当てられる。
【０１０８】
ｆ［１］には、項が電界についてのものであるかまたは磁界についてのものであるかに応じて数値が割り当てられる。ここでは、電界のときを、０とし、磁界のときを、１とする。ｆ［２］には、電界または磁界の座標成分に対応する数値を割り当て、ｘ方向成分であるときは１、ｙ方向成分であるときは２、ｚ方向成分であるときは３とする。
【０１０９】
ｆ［３］には、項の時間に関する微分階数が割り当てられる。例えば、時間の１階微分であるときは１とし、２階微分であるときは２とする。ただし、時間で微分されないときは０とする。同様にして、ｆ［４］には、ｘ方向の微分階数を割り当て、ｆ［５］には、ｙ方向の微分階数を割り当て、ｆ［６］には、ｚ方向の微分階数が割り当てられる。
【０１１０】
ｆ［７］には、この項を誘電率εで除する回数を割り当て、ｆ［８］には、透磁率μで除する回数が割り当てられる。この様にすれば、例えば、電界のｚ方向成分Ｅｚにおけるｘ方向の１階微分かつ時間の１階微分項は、ｆ［０］＝１，ｆ［１］＝０，ｆ［２］＝３，ｆ［３］＝１，ｆ［４］＝１，ｆ［５］＝０，ｆ［６］＝０，ｆ［７］＝０，ｆ［８］＝０として表すことができる。
【０１１１】
以上のようにして項が表され、項の組み合わせにより「式」が表現される。式に含まれる項の数をｓで表し、第１項をｆ１、第２項をｆ２などとする。例えば、∂Ｈｙ／ε∂ｘ−∂Ｈｘ／ε∂ｙという式は、ｓ＝２，ｆ１［０］＝１，ｆ１［１］＝１，ｆ１［２］＝２，ｆ１［３］＝０，ｆ１［４］＝１，ｆ１［５］＝０，ｆ１［６］＝０，ｆ１［７］＝１，ｆ１［８］＝０，ｆ２［０］＝−１，ｆ２［１］＝１，ｆ２［２］＝１，ｆ２［３］＝０，ｆ２［４］＝０，ｆ２［５］＝１，ｆ２［６］＝０，ｆ２［７］＝１，ｆ２［８］＝０として表すことができる。
【０１１２】
そこで、時間微分生成部４は、電界及び磁界の時間の１〜５階微分量と、電界及び磁界の空間微分における時間の１〜４階微分量を、以下の処理手順により算出する。
【０１１３】
時間微分量は、予め項ｆ１により表されているものとする。まず、ｓ＝１，ｍ＝１とし（ステップＳ４１）、項ｆｍを一時的な配列Ｆにコピーする（ステップＳ４２）。ｎ＝０〜８に対して、Ｆ［ｎ］＝ｆｍ［ｎ］である。
【０１１４】
次に、この項ｆｍに時間微分が含まれていれば、マックスウェルの方程式に基づいて展開する。すなわち、Ｆ［３］＞０のときには、ｓ＝ｓ＋１とし、新たな項ｆｓのための記憶領域を確保する（ステップＳ４３，Ｓ４４）。そして、項ｆｍが電界についてのものであるか、磁界についてのものであるかに応じて展開する（ステップＳ４５，Ｓ４６，Ｓ４７）。展開後は、ステップＳ４２に戻る。
【０１１５】
一方、項ｆｍに時間微分が含まれていなければ、ｍ＝ｍ＋１とし（ステップＳ４８）、ｍ＝ｓでなければ、ステップＳ４２に戻る。ｍ＝ｓであれば、展開結果に基づいて時間微分量を算出する（ステップＳ４９，Ｓ５０）。
【０１１６】
図１２のステップＳ５１〜Ｓ５６は、本発明の実施の形態４による電磁界解析装置の動作の一例を示したフローチャートであり、微分演算における電界または磁界の展開手順を表している。微分演算における電界の展開は、以下のようにして行われる。
【０１１７】
まず、符号の演算が行われ、ｆｍ［０］＝Ｆ［０］，ｆｓ［０］＝−Ｆ［０］とする（ステップＳ５１）。
【０１１８】
次に、磁界への変換が行われ、Ｆ［２］＝１のときには、式（５７ａ）に基づいて、ｆｍ［１］＝１，ｆｓ［１］＝１，ｆｍ［２］＝３，ｆｓ［２］＝２とし、Ｆ［２］＝２のときには、式（５７ｂ）に基づいて、ｆｍ［１］＝１，ｆｓ［１］＝１，ｆｍ［２］＝１，ｆｓ［２］＝３とし、Ｆ［２］＝３のときには、式（５７ｃ）に基づいて、ｆｍ［１］＝１，ｆｓ［１］＝１，ｆｍ［２］＝２，ｆｓ［２］＝１とする（ステップＳ５２）。
【０１１９】
そして、時間微分の階数の演算が行われ、ｆｍ［３］＝Ｆ［３］−１，ｆｓ［３］＝Ｆ［３］−１とする（ステップＳ５３）。
【０１２０】
また、空間微分の階数の演算が行われる。Ｆ［２］＝１のときには、式（５７ａ）に基づいて、ｆｍ［４］＝Ｆ［４］，ｆｓ［４］＝Ｆ［４］，ｆｍ［５］＝Ｆ［５］＋１，ｆｓ［５］＝Ｆ［５］，ｆｍ［６］＝Ｆ［６］，ｆｓ［６］＝Ｆ［６］＋１、Ｆ［２］＝２のときには、式（５７ｂ）に基づいて、ｆｍ［４］＝Ｆ［４］，ｆｓ［４］＝Ｆ［４］＋１，ｆｍ［５］＝Ｆ［５］，ｆｓ［５］＝Ｆ［５］，ｆｍ［６］＝Ｆ［６］＋１，ｆｓ［６］＝Ｆ［６］、Ｆ［２］＝３のときには、式（５７ｃ）に基づいて、ｆｍ［４］＝Ｆ［４］＋１，ｆｓ［４］＝Ｆ［４］，ｆｍ［５］＝Ｆ［５］，ｆｓ［５］＝Ｆ［５］＋１，ｆｍ［６］＝Ｆ［６］，ｆｓ［６］＝Ｆ［６］とする（ステップＳ５４）。
【０１２１】
また、ε及びμで除する回数の演算が行われる。ｆｍ［７］＝Ｆ［７］＋１，ｆｓ［７］＝Ｆ［７］＋１，ｆｍ［８］＝Ｆ［８］，ｆｓ［８］＝Ｆ［８］とする（ステップＳ５５，Ｓ５６）。
【０１２２】
一方、磁界の展開も、電界の展開と同様にして行われる。まず、符号の演算が行われ、ｆｍ［０］＝−Ｆ［０］，ｆｓ［０］＝Ｆ［０］とする（ステップＳ５１）。
【０１２３】
次に、電界への変換が行われ、Ｆ［２］＝１のときには、式（５８ａ）に基づいて、ｆｍ［１］＝０，ｆｓ［１］＝０，ｆｍ［２］＝３，ｆｓ［２］＝２、Ｆ［２］＝２のときには、式（５８ｂ）に基づいて、ｆｍ［１］＝０，ｆｓ［１］＝０，ｆｍ［２］＝１，ｆｓ［２］＝３、Ｆ［２］＝３のときには、式（５８ｃ）に基づいて、ｆｍ［１］＝０，ｆｓ［１］＝０，ｆｍ［２］＝２，ｆｓ［２］＝１とする（ステップＳ５２）。
【０１２４】
時間微分及び空間微分の階数の演算は、電界の展開の場合におけるステップＳ５３及びステップＳ５４と同様にして行われる。そして、ε及びμで除する回数の演算が行われる。ｆｍ［７］＝Ｆ［７］，ｆｓ［７］＝Ｆ［７］，ｆｍ［８］＝Ｆ［８］＋１，ｆｓ［８］＝Ｆ［８］＋１とする（ステップＳ５５，Ｓ５６）。
【０１２５】
例えば、電界Ｅｚのｘ方向１階微分かつ時間１階微分項は、ｓ＝１，ｆ１［０］＝１，ｆ１［１］＝０，ｆ１［２］＝３，ｆ１［３］＝１，ｆ１［４］＝１，ｆ１［５］＝０，ｆ１［６］＝０，ｆ１［７］＝０，ｆ１［８］＝０と表されるので、この項ｆ１を上述した手順にしたがって展開すると、ステップＳ４３において、Ｆ［３］＞０だから、ステップＳ４４で新たに項ｆ２が生成される。そして、ステップＳ４５において、Ｆ［１］＝０だから、ステップＳ４６により展開される。
【０１２６】
従って、Ｆ［２］＝３だから、ｓ＝２，ｆ１［０］＝１，ｆ１［１］＝１，ｆ１［２］＝２，ｆ１［３］＝０，ｆ１［４］＝２，ｆ１［５］＝０，ｆ１［６］＝０，ｆ１［７］＝１，ｆ１［８］＝０，ｆ２［０］＝−１，ｆ２［１］＝１，ｆ２［２］＝１，ｆ２［３］＝０，ｆ２［４］＝１，ｆ２［５］＝１，ｆ２［６］＝０，ｆ２［７］＝１，ｆ２［８］＝０となる。この様にして、各項において時間微分がなくなるまで、展開が繰り返され、全ての項において時間微分がなくなると、展開は終了する。
【０１２７】
図１３のステップＳ６１〜Ｓ６９は、本発明の実施の形態４による電磁界解析装置の動作の一例を示したフローチャートであり、微分演算における時間微分量の計算手順を表している。式（５９）により電界のｘ方向成分が表されるとき、例えば、この電界成分のｘ方向１階微分かつｙ方向１階微分は、次式（６５）により表される。
【数４０】
2004239784.tif 000042
【０１２８】
第ｉ，ｊ，ｋ番目の格子点においては、（ｘ，ｙ，ｚ）＝（０，０，０）として、式（６５）から次式（６６）が得られる。
【数４１】
2004239784.tif 000043
【０１２９】
従って、この格子点上における電界成分の微分量は、式（６６）から式（６０ａ）及び（６０ｂ）に基づいて求めることができる。つまり、補間式を微分する方向に偏微分したものの積により、微分量があらわされる。
【０１３０】
そこで、時間微分量の計算は、以下のようにして行われる。補間式の係数は、予め求められているとする。まず、ｍ＝１，ｖ＝０とし（ステップＳ６１）、ｎ＝０〜８に対して、Ｆ［ｎ］＝ｆｍ［ｎ］とする（ステップＳ６２）。
【０１３１】
次に、Ｆ［１］及びＦ［２］に基づいて、対応する電界または磁界成分の補間式の係数を呼び出す（ステップＳ６３）。ここで、ｐ＝Ｆ［０］／（εのＦ［７］乗×μのＦ［８］乗）とする（ステップＳ６４）。
【０１３２】
そして、Ｆ［４］＝０のとき、ｐ＝ｐとし、Ｆ［４］＝１のとき、ｐ＝ｐ×Ｇｘとし、Ｆ［４］＝２のとき、ｐ＝ｐ×２Ｄｘとし、Ｆ［４］＝３のとき、ｐ＝ｐ×６Ｃｘとし、Ｆ［４］＝４のとき、ｐ＝ｐ×２４Ｂｘとし、Ｆ［４］＝５のとき、ｐ＝ｐ×１２０Ａｘとし、Ｆ［４］≧６のとき、ｐ＝０とする（ステップＳ６５）。ここで、電界または磁界のｘ方向１階微分量をＧｘとした。
【０１３３】
また、Ｆ［５］＝０のとき、ｐ＝ｐとし、Ｆ［５］＝１のとき、ｐ＝ｐ×Ｇｙとし、Ｆ［５］＝２のとき、ｐ＝ｐ×２Ｄｙとし、Ｆ［５］＝３のとき、ｐ＝ｐ×６Ｃｙとし、Ｆ［５］＝４のとき、ｐ＝ｐ×２４Ｂｙとし、Ｆ［５］＝５のとき、ｐ＝ｐ×１２０Ａｙとし、Ｆ［５］≧６のとき、ｐ＝０とする（ステップＳ６６）。ここで、電界または磁界のｙ方向１階微分量をＧｙとした。
【０１３４】
また、Ｆ［６］＝０のとき、ｐ＝ｐとし、Ｆ［６］＝１のとき、ｐ＝ｐ×Ｇｚとし、Ｆ［６］＝２のとき、ｐ＝ｐ×２Ｄｚとし、Ｆ［６］＝３のとき、ｐ＝ｐ×６Ｃｚとし、Ｆ［６］＝４のとき、ｐ＝ｐ×２４Ｂｚとし、Ｆ［６］＝５のとき、ｐ＝ｐ×１２０Ａｚとし、Ｆ［６］≧６のとき、ｐ＝０とする（ステップＳ６７）。ここで、電界または磁界のｚ方向１階微分量をＧｚとした。
【０１３５】
ｍ＝ｓでなければ、ｍ＝ｍ＋１，ｖ＝ｖ＋ｐとし（ステップＳ６８，Ｓ６９）、ステップＳ６２〜Ｓ６９を繰り返す。そして、全ての項について微分量が計算されると、計算結果として時間微分量ｖ＋ｐが算出される。
【０１３６】
本実施の形態によれば、時間微分量算出のための微分演算自体も行わせることができる。
【０１３７】
実施の形態５．
図１４は、本発明の実施の形態５による電磁界解析装置における格子点ごとの電界及びその空間微分並びにそれらの時間微分の配置の一例を示した説明図である。本実施の形態では、電界または磁界のいずれか一方が、波動方程式などを解くことによって予め求められている場合に、他方がマックスウェルの方程式から算出される。例として、実施の形態１と同様に平面波の伝搬解析について説明する。より一般的な電磁界の場合も同様にして行われる。
【０１３８】
ｚ軸方向の格子点上に初期値として、電界のｘ方向成分Ｅｘｉ及びそのｚ方向微分量Ｅｘｚｉと、それらの時間微分量Ｅｘｔｉ及びＥｘｚｔｉを配置し、初期値記憶部２に記憶する。電界のｙ方向成分についても同様であるとする。
【０１３９】
そこで、Ｅｘ及びＥｘｔを格子点間において補間する補間式を、次式（６７）及び（６８）によりそれぞれ定義する。
【数４２】
2004239784.tif 000044
【０１４０】
式（６７）及び（６８）において、係数ＡＥ〜ＤＥ及びＡＥｔ〜ＤＥｔは、補間式係数生成部３により求められる。また、これらＥｘ及びＥｘｔの時間微分は、波動方程式などを解くことによって予め求められているものとする。この様な場合に、磁界のｙ方向成分Ｈｙの時間発展が、式（４９）に基づいて算出される。式（４９）における時間微分量は、以下のようにして求められる。
【０１４１】
Ｈｙの時間１階微分は、マックスウェルの方程式（２１ｂ）から次式（６９）により表される。
【数４３】
2004239784.tif 000045
【０１４２】
従って、式（６７）をｚ方向微分することにより、式（６９）から次式（７０）が得られる。
【数４４】
2004239784.tif 000046
【０１４３】
時間の２階微分は、式（６９）を時間微分することにより得られる次式（７１）から求められる。
【数４５】
2004239784.tif 000047
【０１４４】
従って、式（７１）から次式（７２）が得られる。
【数４６】
2004239784.tif 000048
【０１４５】
同様にして、時間の３〜５階微分が求められ、次式（７３）〜（７５）によりそれぞれ表される。
【数４７】
2004239784.tif 000049
【０１４６】
結局、Ｈｙの時間発展は、式（７０），（７２）〜（７５）からｚ＝０として、次式（７６）に基づいて算出される。
【数４８】
2004239784.tif 000050
【０１４７】
本実施の形態によっても、電磁界の伝搬現象を高精度に解析することができる。
【０１４８】
【発明の効果】
以上、説明したように本発明による電磁界解析装置及び電磁界解析方法によれば、空間微分及び時間発展に関する計算精度が向上し、電磁界の空間分布の時間変化をより高精度に解析することができる。
【図面の簡単な説明】
【図１】本発明の実施の形態１による電磁界解析装置の一構成例を示したブロック図である。
【図２】座標空間に配置された格子点ごとの物理量に対して物理量を格子点間において補間する代数曲線を示した説明図である。
【図３】座標空間における格子点ごとの電界及び磁界の配置を示した説明図である。
【図４】座標空間における格子点ごとの電界及び磁界並びにそれらの空間微分の配置の一例を示した説明図である。
【図５】図１の電磁界解析装置における電磁界解析の動作の一例を示したフローチャートである。
【図６】図１の電磁界解析装置における電磁界解析の動作の一例を示したフローチャートである。
【図７】図１の電磁界解析装置による平面矩形波の伝搬現象の解析結果の一例を示した状態図である。
【図８】従来のＦＤＴＤ法による平面矩形波の伝搬現象の解析結果の一例を示した状態図である。
【図９】本発明の実施の形態２による電磁界解析装置における格子点ごとの電界及び磁界並びにそれらの空間微分の配置の一例を示した説明図である。
【図１０】本発明の実施の形態３による電磁界解析装置の動作の一例を示したフローチャートである。
【図１１】本発明の実施の形態４による電磁界解析装置の動作の一例を示したフローチャートである。
【図１２】本発明の実施の形態４による電磁界解析装置の動作の一例を示したフローチャートである。
【図１３】本発明の実施の形態４による電磁界解析装置の動作の一例を示したフローチャートである。
【図１４】本発明の実施の形態５による電磁界解析装置における格子点ごとの電界及びその空間微分並びにそれらの時間微分の配置の一例を示した説明図である。
【符号の説明】
１  電磁界解析装置、２  初期値記憶部、３  補間式係数生成部、
４  時間微分生成部、５  時間発展算出部、６  ディスプレイ、７  補助記憶部
【図１】
2004239784.tif 000051
【図２】
2004239784.tif 000052
【図３】
2004239784.tif 000053
【図４】
2004239784.tif 000054
【図５】
2004239784.tif 000055
【図６】
2004239784.tif 000056
【図７】
2004239784.tif 000057
【図８】
2004239784.tif 000058
【図９】
2004239784.tif 000059
【図１０】
2004239784.tif 000060
【図１１】
2004239784.tif 000061
【図１２】
2004239784.tif 000062
【図１３】
2004239784.tif 000063
【図１４】
2004239784.tif 000064