JP2008250947A - 波動場解析方法および装置ならびに波動場解析方法をコンピュータに実行させるプログラムを記録したコンピュータに読み取り可能な記録媒体 - Google Patents

Abstract

【課題】複数の、類似した複数の計算モデルを繰返し計算する場合に、計算量を低減し、高速化することができる波動場解析方法、装置および波動場解析方法をコンピュータに実行させるプログラムを記録したコンピュータに読み取り可能な記録媒体を提供する。
【解決手段】解析領域内の波動の伝播状態が異なる複数の計算モデルにおいて、初期状態から波動の伝播状態が共通である共通部を抽出し、抽出された共通部について、波動場の時間変化の計算を初期状態から開始し、途中の計算段階において計算を停止し、計算された波動場の状態を記憶手段に保存し、記憶手段に保存された波動場の状態を読み出し、途中の計算段階の波動場の状態から、各モデルに対する波動場の時間変化の計算を連続して行うことにより、上記課題を解決する。
【選択図】図５

Description

本発明は、波動場解析方法および装置ならびに波動場解析方法をコンピュータに実行させるプログラムを記録したコンピュータに読み取り可能な記録媒体に関し、詳しくは、ＦＤＴＤ法（Finite Difference Time Domain method（時間領域差分法））を用いて、電磁波などの波動の伝播、すなわち電磁場（電磁界）などの波動場の時間的変化を数値解析する波動場解析方法および波動場解析装置ならびに波動場解析方法をコンピュータに実行させるプログラムを記録したコンピュータ読み取り可能な記録媒体に関する。

近年、光ディスク、レンズ、特にマイクロレンズなどの種々の光学部品、光コネクタ、光導波路、光ファイバなどの光通信用部品、アンテナ、レーザ装置、顕微鏡、半導体集積回路（ＩＣ）、圧電素子、構造物等々の設計のために、流体内や固体内の波動解析シミュレーション法として、例えば、光やＸ線などを含む電磁波、音波や、弾性波なども含め、波動の伝播の時間的な変化を解析するシミュレーション法として、ＦＤＴＤ法が利用されてきている。

このＦＤＴＤ法は、特に、マックスウェル（Maxwell）方程式に従う電磁波の伝播の時間的な変化の解析ツールとして利用されている。このＦＤＴＤ法を用いて、対象とする伝播系、すなわち、解析領域内において、マックスウェル方程式をコンピュータで数値的に解くことにより、この解析領域内において、例えば空気中を電磁場が伝播し、対象とする解析領域内に存在する物質表面で反射したり、屈折して透過したりする伝播状態を理論的に予測することができる。
このような電磁波の伝播状態の論理的な予測は、上述の光学部品、光通信用部品、アンテナを含む種々の部品等の設計に極めて有用な情報や必要な情報を与える。

このようなＦＤＴＤ法を用いた電磁場解析装置および方法が、特許文献１に提案されている。この特許文献１には、解析領域の全領域（全空間）の電磁場がゼロ「０」の状態を初期状態として、解析領域の一部分に、入射光を発生する励震部を設定し、励震部から連続波、パルス波などを発生させ、ＦＤＴＤ法によって、解析領域内の電磁場の時間変化を計算することが開示されている。また、特許文献１では、入射波に関する情報を入力し、入力された情報に基づいて入射波を計算することにより入射波の誤差を正確に再現し、再現された誤差を含んだ値を計算された入射波とし、この計算された入射波に基づいて電磁場を計算することによりＦＤＴＤ法による誤差を相殺し、ＦＤＴＤ法による誤差のない計算結果を出力することができることを開示している。

この特許文献１に開示されるようなＦＤＴＤ法を用いた従来の電磁場解析方法のフローチャートを図９に示す。
まず、ステップＳ２００において、異なる複数の計算モデル（図示例では、３つの計算モデル１，２，３）に対して、それぞれ屈折率の空間分布を定義する。
また、ステップＳ２０２において、解析領域、入射波（入力電磁波）およびその他、電磁場計算に必要となる条件を設定するために計算条件を入力する。
ここで、屈折率の空間分布定義ステップＳ２００と、計算条件の入力ステップＳ２０２とは、いずれを先に行っても良いし、同時に行っても良い。
こうして、まず、電磁場計算を開始する前に、各計算モデルについて、電磁場計算の対象となる解析領域の屈折率の空間分布の定義および必要な計算条件の設定がなされる。

その後、ステップＳ２０４において、ｉ＝１として、計算モデル１について、解析領域の屈折率の空間分布の定義および必要な計算条件の設定に基づいて、ＦＤＴＤ法計算プログラムを用いて、電磁場が存在しない状態を初期状態とし、解析空間上に定めた励震部からの電磁波を入射させ、マックスウェル方程式をコンピュータで数値的に解く電磁場の計算を行い、電磁場の時間的、空間的変化を計算する。
続いて、ステップＳ２０６において、計算後処理として、こうして計算された電磁場の時間、空間分布に基づいて、解析目的に応じた物理量を計算する。また、計算結果を、表示デバイスのディスプレー上にグラフィック表示したり、プリンタにハードコピー出力する。
こうして、計算モデル１についての電磁場計算が終了する。

次に、ステップＳ２０８において、全計算モデルの計算が終了したか否かを判断する。図示例では、ｉ＝３（または３以上）であるか否かを判断する。全計算モデルの計算が終了していれば（図示例では、ｉが３以上であれば）、Ｙｅｓとして、ステップＳ２１０（ＥＮＤ）に抜け出し、電磁場解析を終了する。一方、全計算モデルの計算が終了していなければ（図示例では、ｉが３未満であれば）、Ｎｏ：ｉ＝ｉ＋１として、ステップＳ２０４に戻る。

続いて、ステップＳ２０４において、計算モデル２（またはｉ＋１）について、解析領域の屈折率の空間分布の定義および必要な計算条件の設定に基づいて、ＦＤＴＤ法計算プログラムを用いて、最初の状態（初期状態）から電磁場の時間的、空間的変化を計算し、次に、ステップＳ２０６において、計算後処理として、解析目的に応じた物理量の計算や、計算結果の出力を行い、次に、ステップＳ２０８において、全計算モデルの計算が終了したか否かを判断し、終了していなければステップＳ２０４に戻り、全計算モデルの計算が終了するまで、上記ステップＳ２０４〜Ｓ２０８を繰り返し、終了していればステップＳ２１０に抜け出して、電磁場解析を終了する。
このように、従来の電磁場解析技術では、単純に複数の計算モデルに対する計算を、初期状態から、逐次実施している。

特開２０００−２２７４５０号公報

ところで、特許文献１に開示されたような技術では、上述したように、同じ解析領域において類似した電磁波の伝播モデルについて、ＦＤＴＤ法によって電磁波の伝播の時間変化を計算する場合においても、いずれの計算モデルにおいても、初期状態から計算しなおす必要があり、同一の計算を繰り返す必要があるという問題がある。
例えば、マイクロレンズ等の開発や設計のために、レンズ構造の一部、例えば形成被膜の膜厚や、マイクロレンズの曲率などを少しずつ変化させて、光の透過状態などの光の伝播状態の解析を行い、最適なレンズ構造やマイクロレンズの最適な曲率などを探す場合においても、すなわち、計算モデルが如何に類似していたとしても、いずれの類似計算モデルにおいても、最初の初期状態から、計算しなおす必要があり、同一の計算を繰返す必要があるという問題がある。

本発明の課題は、上記従来技術の問題点を解消し、複数の計算モデル、特に、類似した複数の計算モデルを繰返し計算する場合に、計算量を低減し、高速化することができる波動場解析方法および装置ならびに波動場解析方法をコンピュータに実行させるプログラムを記録したコンピュータに読み取り可能な記録媒体を提供することにある。

上記課題を解決するために、本発明の第１の態様は、解析領域内を波動が伝播する伝播状態を解析する波動場解析方法であって、同一の前記解析領域内の前記波動の伝播状態が異なる複数の計算モデルにおいて、初期状態から前記波動の伝播状態が共通である共通部を抽出し、抽出された前記共通部について、前記波動が伝播する前記解析領域内の波動場の時間変化の計算を初期状態から開始し、計算開始後、前記共通部において、前記解析領域内の波動場の時間変化を連続して計算し、前記共通部内の所定の途中の計算段階において、前記計算を停止して、計算された前記波動場の状態を記憶手段に保存し、前記複数の計算モデルの各モデルに対する、前記解析領域内の前記波動場の時間変化の計算を開始するに際し、前記記憶手段に保存された前記共通部内の前記所定の計算段階の前記波動場の状態を読み出し、前記記憶手段から読み出された前記共通部の前記所定の途中の計算段階の前記波動場の状態から、各モデルに対する、前記解析領域内の前記波動場の時間変化の計算を連続して行うことを特徴とする波動場解析方法を提供するものである。

ここで、前記共通部内の所定の途中の計算段階の前記波動場の状態は、前記共通部の最後の計算段階の前記波動場の状態であるのが好ましい。
また、前記波動は、前記解析領域内に設定された励震部（波源）において発生され、前記共通部は、前記励震部（波源）近傍の前記波動を伝播する媒質が共通である領域であるのが好ましい。
また、前記波動は、電磁波であり、前記波動場は、電磁場であるのが好ましい。
また、前記電磁波は、前記解析領域内に設定された励震部（波源）において発生され、前記共通部の抽出は、前記解析領域内に前記初期状態として初期の電磁場が存在せず、前記励震部（波源）近傍の屈折率分布が共通である領域を抽出することであり、前記共通部内の所定の途中の計算段階は、電磁波が前記共通部内を伝搬する時間に設定されるのが好ましい。
また、前記複数の計算モデルが、前記電磁場の計算結果に基いて修正されるのが好ましい。

また、上記課題を解決するために、本発明の第２の態様は、解析領域内に波動が伝播する伝播状態を解析する波動場解析装置であって、前記波動が伝播する前記解析領域内の波動場の時間変化を初期状態から連続して計算する波動場計算手段と、同一の前記解析領域内の前記波動の伝播状態が異なる複数の計算モデルにおいて、前記初期状態から前記波動の伝播状態が共通である共通部を抽出する共通部抽出手段と、前記共通部抽出手段によって抽出された前記共通部について、前記波動場計算手段によって前記波動場の時間変化を前記初期状態から連続して計算することにより得られた、前記共通部内の所定の途中の計算段階の前記波動場の状態を保存する記憶手段と、前記記憶手段に保存された前記共通部内の前記所定の途中の計算段階の前記波動場の状態を読み出す読出手段とを有し、前記読出手段によって前記記憶手段から読み出された前記共通部の前記所定の途中の計算段階の前記波動場の状態から、前記波動場計算手段によって、前記複数の計算モデルの各モデルに対する、前記解析領域内の前記波動場の時間変化の計算を開始することを特徴とする波動場解析装置を提供するものである。

ここで、前記共通部内の所定の途中の計算段階の前記波動場の状態は、前記共通部の最後の計算段階の前記波動場の状態であるのが好ましい。
また、さらに、前記解析領域内に設定され、前記波動を発生させる励震部（波源）を有し、
前記共通部抽出手段は、前記共通部として、前記励震部（波源）近傍の前記波動を伝播する媒質が共通である領域を抽出するのが好ましい。
また、前記波動は、電磁波であり、前記波動場は、電磁場であるのが好ましい。
また、さらに、前記解析領域内に設定され、前記電磁波を発生させる励震部（波源）を有し、
前記共通部抽出手段は、前記解析領域内に前記初期状態として、初期の電磁場が存在せず、前記励震部（波源）近傍の屈折率分布が共通である領域を抽出するものであり、
前記共通部内の所定の途中の計算段階は、電磁波が前記共通部内を伝搬する時間に設定されるのが好ましい。
また、前記複数の計算モデルが、前記電磁場の計算結果に基いて修正されるのが好ましい。

また、上記課題を解決するために、本発明の第３の態様は、上記第１の態様の波動場解析方法をコンピュータに実行させるプログラムを記録したことを特徴とするコンピュータに読み取り可能な記録媒体を提供するものである。

本発明によれば、上記構成により、同一の解析領域内の電磁波などの波動の伝播状態が異なる複数の計算モデルであって、初期状態から波動の伝播状態が共通である領域（共通領域）がある類似した計算モデルを計算する場合に、まず、最初に、初期状態から、共通領域を伝播するまでの時間後の、途中段階の電磁場などの波動場を保存しておき、異なる複数の計算モデルの各計算モデルを計算する場合には、途中段階の波動場から計算を開始することにより、複数の計算モデルに対し、途中までの計算を省略できるので、計算時間を短縮することができる。

以下に、本発明に係る波動場解析方法および装置ならびに波動場解析方法をコンピュータに実行させるプログラムを記録したコンピュータに読み取り可能な記録媒体を、添付の図面に示す好適実施形態に基づいて詳細に説明する。

本発明の波動場解析方法および装置は、同一の解析領域内の電磁波等の波動の伝播状態が異なる複数の計算モデルをマックスウェル方程式などの波動方程式を用いて計算する場合、特に、初期（開始）状態から波動の伝播状態が共通であるモデル領域（共通領域）、例えば、励震部（波源）近傍に共通のモデル領域がある類似した計算モデルを計算する場合に、まず最初に、初期状態から、共通領域を伝播するまでの時間後の、途中段階の電磁場などの波動場を保存しておき、異なる複数の計算モデルの各計算モデルを計算する場合には、途中段階の波動場から計算を開始することにより、類似した計算モデルを繰返し計算する場合に、計算量を低減して、高速化するものである。

以下では、波動の代表例として、光や電子線やＸ線などを含む電磁波を用い、波動場および波動方程式の代表例として、電磁場およびマックスウェル方程式を用い、ＦＤＴＤ法によってマックスウェル方程式をコンピュータで数値的に解き、対象とする解析領域の電磁波の空間および時間分布を計算し、電磁場の時間的な変化を求めて電磁場の解析を行う例を代表例として説明するが、本発明はこれに限定されるわけではない。本発明は、波動方程式をコンピュータで数値的に解き、対象とする解析領域の波動の空間および時間分布を計算し、波動場の時間的な変化を求めて波動場の解析を行うものであれば、例えば、音波や、弾性波などの電磁波以外の波動にも適用可能であることはもちろんである。

すなわち、本発明は、電磁波の空間、時間分布を、マックスウェル方程式に基き、屈折率（誘電率）の空間分布を境界条件として、コンピュータで数値的に解くためのＦＤＴＤ法プログラムを用いたシミュレーション方法および手段を提供するもので、主に、光ディスク、レンズ、特に、マイクロレンズ等の種々の光学部品、光コネクタ、光導波路、光ファイバ等の光通信用部品、アンテナ、レーザ装置、顕微鏡、半導体集積回路（ＩＣ）、圧電素子、構造物等々の解析、微小領域の光学現象解析などに用いられ、開始状態の電磁場を０とし、外部から入射された電磁波の時間変化を計算する手法を提供するものである。

まず、本発明の一実施形態に係る波動場解析方法を実施する波動場解析装置が適用される電磁場解析装置の機能的構成について説明する。
図１は、本発明の波動場解析装置が適用される電磁場解析装置の機能的構成の一実施例を示すブロック図である。
図１に示す本発明の電磁場解析装置１０は、パソコンやワークステーションなどの計算機（コンピュータ）環境で動作するプログラムを含み、コンピュータによって構成することができるものである。

図１に示すように、電磁場解析装置１０は、複数の計算モデルに対する屈折率の空間分布を定義する屈折率分布定義部１２と、解析領域を含む電磁場計算に必要となる条件を設定するために計算条件を入力する計算条件入力部１４と、屈折率分布定義部１２によって定義された計算モデルの屈折率分布に基づいて、解析領域内の電磁波の伝播状態が共通となる共通部を抽出する共通部抽出部１６と、上述の計算モデルの屈折率分布および上述の計算条件に基づいて、ＦＤＴＤ法を用いてマックスウエル方程式を数値的に解く電磁場計算部１８と、電磁場計算部１８によって計算された途中段階の電磁場（空間分布）などを記憶保存する記憶部２０と、途中段階の電磁場のデータの記憶部２０への書き込みおよび記憶部２０からの途中段階の電磁場データの読み出しを制御して行う書込／読出制御部２２と、電磁場計算部１８によって、記憶部２０から書込／読出制御部２２によって読み出された途中段階の電磁場から開始され、残りの解析領域の末端で終了した全空間を伝播した電磁場の計算結果を出力する計算結果出力部２４とを有する。

屈折率分布定義部１２は、解析対象となる解析領域の形状、材料に基づいて、計算モデルとなる屈折率の空間分布を定義するものである。屈折率の空間分布は、異なる複数の計算モデルに対してそれぞれ定義される。本発明では、複数の計算モデルは類似しているものを対象としているので、複数の計算モデルの各々に対して定義された解析領域内の屈折率の空間分布には、屈折率が同一である空間領域が存在する。
計算条件入力部１４は、計算範囲、空間離散化方法、境界条件などの解析領域の入力や設定、空間分布、時間分布、偏光条件などの入力電磁波（入射波）の入力や設定、およびその他、時間離散化条件などの計算に必要となる条件の入力や設定を行うものである。

これらの屈折率分布定義部１２および計算条件入力部１４は、電磁場計算部１８によって解析領域内の電磁場計算を開始する前に、複数の計算モデルの各々の計算モデルについて、電磁場計算の対象となる解析領域の屈折率の空間分布の定義および必要な計算条件の入力・設定を行う。
なお、屈折率分布定義部１２による屈折率の空間分布定義と、計算条件入力部１４による計算条件の入力および／または設定とは、互いに対応させながら同時に行っても良いし、いずれかを先に行っても良いし、また、独立して行っても良い。

共通部抽出部１６は、本発明の特徴とする部分であって、屈折率分布定義部１２によって定義された計算モデルの屈折率分布、またはこの屈折率分布および計算条件入力部１４によって設定された計算条件に基づいて、解析領域内の電磁波の伝播状態が共通となる領域、すなわち、屈折率分布定義部１２で定義された解析領域内の屈折率の空間分布に存在する、屈折率が同一である空間領域を、共通部（共通領域）として抽出するものである。
共通部抽出部１６による共通部の抽出は、電磁場計算部１８による解析領域内の電磁場計算の開始前に行われる。なお、共通部の抽出の詳細については、後述する。

電磁場計算部１８は、屈折率分布定義部１２によって定義された計算モデルの屈折率分布および計算条件入力部１４によって設定された計算条件に基づいて、ＦＤＴＤ法計算プログラムによって、マックスウエル方程式をパソコンやワークステーションなどのコンピュータで数値的に解くものであり、解析領域（解析空間）内の電磁場の時間的、空間的変化を計算する。
ここで、電磁場計算部１８が実施するＦＤＴＤ法計算プログラムは、電磁場の物理現象を表すマックスウエル方程式を、コンピュータで数値的に解くための１つの手法であり、電磁場が存在しない状態を初期状態とし、解析空間上に定めた励震部（波源）からの電磁波を入射させ、電磁場の時間変化を逐次更新しながら計算するシミュレーション手法である。（なお、ＦＤＴＤ法やＦＤＴＤ法計算プログラムについては、「ＦＤＴＤ法による電磁界およびアンテナ解析」宇野 亨、コロナ社発行の文献を参照することができる。）

ＦＤＴＤ法やＦＤＴＤ法計算プログラムについては、上記文献を始めとして、種々の従来公知の文献に開示された技術を用いれば良いので、詳細な説明は省略するが、以下に、本発明において適用されるＦＤＴＤ法の一例について簡単に説明する。
ＦＤＴＤ法は、下記式（１）および（２）に示すマックスウェル方程式を、時間、空間方向に差分化した電磁場の数値計算手法である。
電磁場計算では、下記式（１）および（２）に示すマックスウェル方程式を、時間および空間を離散化して差分方程式、すなわち、下記式（３）および（４）に示す電磁場計算差分式(一例として電磁場の１成分のみが示されている)に変換する。

本発明に適用されるＦＤＴＤ法においては、図２に示すように、入射光（電場Ｅ、磁場Ｈの電磁波）が入射する解析空間（領域）をメッシュ状に分割して離散化し、解析空間の形状や材料に応じて、各格子（微小セル）点に、屈折率ｎ（誘電率ε）を割り当てる。図２において、解析空間から１つの格子（単位セル）を抜き出して示すように、Ｙｅｅセルと呼ばれる１つの格子上の点に、電場および磁場のｘ，ｙ，ｚ成分、すなわち、Ｅｘ，Ｅｙ，Ｅｚ，Ｈｘ，Ｈｙ，Ｈｚを定義する。

ここで、上記式（３）から明らかなように、磁場（磁界）の変化（差分）は、周りの電場（電界）の和および差の定数倍、すなわち、図３に示されるように、周りの電場が形成する渦の定数倍で表され、上記式（４）から明らかなように、電場の変化（差分）は、周りの磁場の和および差の定数倍、すなわち、図３に示されるように、周りの磁場が形成する渦の定数倍で表される。
このため、ＦＤＴＤ法では、各Ｙｅｅセルにおいて、ある時間の電場から磁場を上記式（３）で計算し、これで求まった磁場を上記式（４）に代入して、次の時間の電場を計算し、さらに、これで求まった電場を再び上記式（３）に代入して、磁場を計算する、あるいは、ある時間の磁場から電場を上記式（４）で計算し、これで求まった電場を上記式（３）に代入して次の時間の磁場を計算し、さらに、これで求まった磁場を再び上記式（４）に代入して、電場を計算する。

ＦＤＴＤ法においては、解析空間の全Ｙｅｅセルについて計算を繰り返し行うことにより、電磁場の時間変化を計算することができる。まず、ある時間の、解析空間における電場の空間分布から、半ステップ後の磁場の空間分布を計算する。引続き、この解析空間における磁場分布から、半ステップ後の電場を計算することで、1ステップ後の電場に更新される。このように、電場、磁場を交互に更新しながら、電磁場の時間変化を計算する。
したがって、ＦＤＴＤ法では、このように、電場と磁場とを交互に計算する計算手順（リープフロッグ（蛙飛び）アルゴリズム）を繰り返すことにより、電磁場の計算を時系列に沿って行うことができる。

以上から、本発明では、ＦＤＴＤ法において、解析領域内に励震源（波源）を定め、入射波を与え、解析領域内の全空間の電磁場を０とした状態を初期（開始）状態として、上述した式（３）および（４）による計算を繰り返し行うことにより、入射波が、計算モデル内を伝搬する様子を計算することができる。
また、本発明では、入射波として、連続波を与え、変化が収束するまで計算を続ける、もしくは、パルス波を与えることもできる。また、本発明において、境界としては、吸収境界条件、周期境界条件などを設定することができる。
本発明においては、上述したように、ＦＤＴＤ法やＦＤＴＤ法計算プログラムを適用することができる。

本発明において、電磁場計算部１８は、まず、共通部抽出部１６によって抽出された、解析領域内の屈折率が同一である共通部について、上述したＦＤＴＤ法を用いて、電磁場が存在しない状態を初期状態とし、解析空間（共通部）上に定めた励震部（波源）からの電磁波を入射させて、電磁場計算を開始し、共通部（共通領域）内の電磁波の伝播に伴う共通部内の電磁場の時間的および空間的変化を計算するものである。
続いて、電磁場計算部１８は、共通部を伝播するまでの時間後の、途中段階の電磁場の計算結果を、書込／読出制御部２２を介して記憶部２０に書き込み、記憶保存させる。

次に、電磁場計算部１８は、いずれの計算モデルにおいても、途中段階の電磁場の計算結果を、記憶部２０から書込／読出制御部２２を介して読み出し、途中段階から残りの解析領域内の電磁波の伝播に伴う電磁場の時間的および空間的変化を、上述したＦＤＴＤ法を用いて、計算するものである。
なお、電磁場計算部１８が記憶部２０に記憶保存させる途中段階の電磁場の計算結果は、共通部であれば、何処まで電磁波が伝播した状態の電磁場の計算結果であっても良いが、共通部の終端まで電磁波が伝播した状態の電磁場の計算結果、すなわち、共通部の最後の計算段階の電磁場の計算結果であるのが好ましい。その理由は、残りの解析領域内の電磁場の計算の量をより少なくすることができ、計算をより高速化できるからである。

記憶部２０は、電磁場計算部１８によって計算された、初期状態から、共通部抽出部１６によって抽出された共通部を伝播するまでの時間後の、途中段階の電磁場などを記憶保存するものである。
書込／読出制御部２２は、電磁場計算部１８によって計算された途中段階の電磁場を記憶部２０に記憶保存させるために途中段階の電磁場のデータを記憶部２０に書き込むと共に、その書き込み制御を行い、かつ、記憶部２０に記憶保存された途中段階の電磁場のデータを記憶部２０から読み出すと共に、その読み出し制御を行うものである。

計算結果出力部２４は、電磁場計算部１８による１つの計算モデルの電磁場計算終了後の後処理として、電磁場計算部１８によって、記憶部２０から書込／読出制御部２２によって読み出された途中段階の電磁場から開始され、残りの解析領域の末端で終了した全空間を伝播した電磁場の時間および空間分布の計算結果に基づいて、反射率、透過率、吸収率および遠方放射分布などの解析目的に応じた物理量を計算する。また、計算結果出力部２４は、こうして計算された物理量や、電磁場の時間および空間分布などの計算結果を、表示デバイスのディスプレー上にグラフィック表示する機能などを備えていても良いし、プリンタでハードコピーを出力する機能などを備えていても良い。

次に、図１に示す電磁場解析装置のハードウエア構成について説明する。
図４は、図１に示す電磁場解析装置のハードウエア構成の一実施例を示すブロック図である。
図４に示す電磁場解析装置３０は、図１に示す電磁場解析装置１０の機能的構成をハードウエア構成として持つもので、上述したように、図１に示す本発明の電磁場解析装置１０を、パソコンやワークステーションなどのコンピュータなどで構成した場合の一実施例である。

図４に示す電磁場解析装置３０は、解析領域や、屈折率や、計算条件などのデータの入力を行う入力装置３２と、本発明の電磁場解析方法をコンピュータなどで実行するプログラム、ＦＤＴＤ法計算プログラムなどのプログラムや、解析領域、屈折率、計算条件などのデータを記憶した着脱可能なデータ記録媒体からデータを読み取り、またデータ記録媒体にデータを書き込む読取／書込装置３４と、ＦＤＴＤ法計算プログラムを実行して電磁場の時間的および空間的な変化を計算し、かつ、共通部を抽出して設定するなどの演算等の処理を行う、また、途中段階での電磁場の時間的および空間的な分布の書き込み／読み出しを初めとして全体を制御する演算処理装置（例えば、ＣＰＵ）３６と、本発明の特徴とする共通部である途中段階の電磁場の時間および空間分布、ＦＤＴＤ法計算プログラムや、本発明の電磁場解析方法の実行プログラムなどのプログラムや、必要な計算条件などのデータを記憶するメモリ３８と、計算結果等を表示する表示装置４０と、計算結果等を記録する記録装置（プリンタ）４２と、これらの各装置を接続するためのバス４４とを有する。

ここで、入力装置３２は、コンピュータ等のキーボードやマウスなどを含み、また、読取／書込装置３４は、各種のデータ記録媒体（メディア）の読み書きをそれぞれ行うメディアドライブなどを含む。また、入力装置３２および読取／書込装置３４は、図１に示す電磁場解析装置１０の複数の計算モデルに対する屈折率の空間分布を定義する屈折率分布定義部１２および解析領域を含む電磁場計算に必要となる条件を設定するために計算条件を入力する計算条件入力部１４としての機能を含む。
また、演算処理装置３６は、例えば、コンピュータ等のＣＰＵ（中央処理ユニット）などで構成され、解析領域内の電磁波の伝播状態が共通となる共通部を抽出する共通部抽出部１６およびＦＤＴＤ法プログラムを実行してマックスウエル方程式を数値的に解く電磁場計算部１８、ならびに途中段階の電磁場のデータの記憶部２０への書き込および記憶部２０からの途中段階の電磁場データの読み出しを制御して行う書込／読出制御部２２としての機能を含む。

メモリ３８は、コンピュータ等のハードディスクや、各種の内蔵または外付けのメモリを含み、途中段階の電磁場（時間、空間分布）などを記憶保存する記憶部２０としての機能を含む。
また、表示装置４０は、計算結果などを表示画面に表示するためのもので、コンピュータなどの液晶表示デバイスやＣＲＴなどのディスプレーやモニタを含み、また、記録装置４２は、計算結果などを紙などに記録したハードコピーを出力するためのもので、いわゆるプリンタを含む。表示装置４０および記録装置（プリンタ）４２は、途中段階の電磁場から開始され、残りの解析領域の末端で終了した全空間を伝播した電磁場の計算結果を出力する計算結果出力部２４としての機能を持つ。

本発明の波動場解析装置が適用される電磁場解析装置は、基本的に以上のような機能的構成およびハードウエア構成を有するものであるが、以下に、上述した電磁場解析装置の作用および電磁場解析の実施例を説明することにより、本発明の波動場解析方法が適用される電磁場解析方法を説明する。
図５は、本発明の波動場解析方法が適用される電磁場解析方法の一例を示すフローチャートである。なお、図５に示す電磁場解析方法は、図１に示す電磁場解析装置によって実施されるものである。
なお、以下では、空気中に存在する球形粒子の散乱の、粒子サイズ依存性を解析する実施例を代表例として説明を行なうが、本発明はこれに限定されないことはいうまでもないことである。

まず、図５に示す電磁場解析方法のステップＳ１００において、図１に示す電磁場解析装置１０の屈折率分布定義部１２によって、解析領域（空間）に対して、それぞれ異なる屈折率の空間分布が定義され、異なる複数の計算モデル（図示例では、３つの計算モデル１，２，３）が定義される。
例えば、図６（ａ），（ｂ）および（ｃ）に示すように、空気（屈折率ｎ＝１．０）中の立方体形状（１０μｍ×１０μｍ×１０μｍ）の３次元空間からなる解析領域５０（５０ａ，５０ｂ，５０ｃ）が定義され、その解析領域内の中心に、水滴（ｎ＝１．３３）からなる粒子サイズの異なる（例えば、直径１μｍ，２μｍ，３μｍの）球形粒子５２（５２ａ，５２ｂ，５２ｃ）が位置する異なる３つの計算モデル１，２および３が定義される。すなわち、ここでは、異なる３つの計算モデル１，２および３に対して、解析領域５０内の中心の直径１μｍ，２μｍ，３μｍの球形水滴粒子５２（５２ａ，５２ｂ，５２ｃ）部分の屈折率ｎが１．３３で、この球形水滴粒子の部分を除いた残りの空気の解析領域５０内の屈折率ｎが１．０となる屈折率分布が定義される。

また、ステップＳ１０２において、計算条件入力部１４によって、解析領域５０の計算範囲として１０μｍ×１０μｍ×１０μｍの３次元空間、空間離散化方法として５０ｎｍの空間グリッド、したがって、Ｙｅｅセルとして５０ｎｍ×５０ｎｍ×５０ｎｍの単位セル（格子）が設定され、入射光（電磁波）条件として、解析領域５０の一端面（図６（ａ），（ｂ）および（ｃ）中、解析領域５０ａ，５０ｂ，５０ｃの上側端面）が励震部５４とし、励震部５４から波長５５０ｎｍの平面波が入射されることが設定され、また、その他の計算条件として、解析領域境界には、ＰＭＬ（完全一致層：perfectly matched layer）吸収境界条件が設定される。こうして、ステップＳ１０２においては、必要な計算条件が入力される。
ここで、ステップＳ１００とステップＳ１０２とは、いずれを先に行っても良いし、同時に行っても良い。これらのステップＳ１００およびＳ１０２によって、電磁場計算を開始する前に、まず、各計算モデルについて、電磁場計算の対象となる解析領域の屈折率の空間分布の定義および必要な計算条件の設定がなされる。

次に、ステップＳ１０４において、共通部抽出部１６によって、計算モデル１，２，３の計算メッシュ上の屈折率を励震部（波源）５４側から比較し、屈折率が同一である共通部を抽出する。例えば、図７（ａ），（ｂ）および（ｃ）に示すように、解析領域５０の励震部５４の面側から、空間グリッドの1層ずつ、屈折率の同一性がなくなるまで、計算モデル１，２，３の各層の屈折率の同一性を判断していき、屈折率が同一である共通部（領域）５６を見つける。
図６および図７に示す例では、計算モデル１，２，３の球形水滴粒子５２ａ，５２ｂおよび５２ｃのサイズが、それぞれ直径１μｍ、２μｍおよび３μｍであり、解析領域５０が、１０μｍ×１０μｍ×１０μｍの３次元空間であるので、励震部５４側の上端面から直径３μｍの球形水滴粒子５２ｃに至らない位置、３．５μｍ（＝Ｌ）の面までの空気中の空間領域の屈折率が同一（ｎ＝１．０）である空間的な共通部５６である。

したがって、共通部５６は、離散化空間グリッドが５０ｎｍであることから、３．５μｍ／５０ｎｍ＝７０層となる。
その結果、このステップＳ１０４においては、７０層分の光路長を電磁波が伝搬する時間Ｔは、ｃを光速（ｃ＝３．００×１０^８ｍ／ｓ）とする時、Ｔ＝Ｌ／ｃ（＝３．５μｍ／ｃ＝３．５×１０^−６ｍ／３．００×１０^８ｍ／ｓ＝１１．７×１０^−１５ｍ／ｓ＝１１．７ｆｓ）となるので、この伝播時間Ｔ（＝１１．７ｆｓ）を、ＦＤＴＤ法計算の時間ステップΔｔ（＝０．０７６４４ｆｓ／ｓｔｅｐ）で割ったＴ/Δｔを、共通部５６の伝搬ステップ数Ｎ（＝１５３ｓｔｅｐ）として設定し、記憶する。

次に、ステップＳ１０６において、ｉ＝１として、図６（ａ）に示す計算モデル１について、電磁場計算部１８によって、解析領域５０の屈折率の空間分布の定義、必要な計算条件および共通部５６の設定に基づいて、ＦＤＴＤ法計算プログラムを用いて、電磁場が存在しない状態を初期状態とし、解析領域の上端面に定めた励震部５４からの入射光を入射させ、マックスウェル方程式をコンピュータで数値的に解き、共通部５６の電磁場の計算を行い、共通部５６の電磁場の時間的、空間的変化を計算する。
なお、電磁場計算部１８による電磁場の計算では、解析領域５０の全空間の電磁場が０の状態を開始状態として、計算を進める。電磁場計算部１８では、励震部５４で与えた入射光（電磁波）が、徐々に解析空間５０内を伝播して行くので、球形粒子５２に当たっても、球形粒子５２との散乱を計算できる。電磁場計算部１８では、一定の計算ステップ、もしくは、計算の収束状態の判定条件など終了条件を満たした段階で計算を終了させる。

このステップＳ１０６では、計算開始から、共通部５６、すなわち、励震部５４付近の共通の構造部分を伝播している間は、共通部５６の中での電磁場は、計算モデルによらず同じである。
そのため、ステップＳ１０６においては、電磁場の計算を、共通部５６の伝搬ステップ数Ｎ（＝１５３ｓｔｅｐ）だけ行い、共通部５６を伝播した電磁場の空間分布、すなわち共通部５６の終端面（図中の下端面）についての、共通部５６の最後の計算段階（最後の時間ステップ（１５３ｓｔｅｐ目））の電磁場の空間分布を求める。

次に、ステップＳ１０８において、図６（ａ）に示す計算モデル１についてステップＳ１０６で計算された共通部５６の最後の計算段階の時間ステップ（Ｎ＝１５３ｓｔｅｐ目）の電磁場の空間分布を、読出／書込制御部２２によって記憶部２０に書き込んで、記憶保存する。
換言すれば、このステップＳ１０６では、上述したように、計算モデル1に対し、通常通り、全空間の電磁場が０の状態から計算が開始され、ステップＳ１０４で設定した、共通部５６の伝搬ステップ数Ｎ（＝１５３ｓｔｅｐ）の計算段階で計算を中断し、ステップＳ１０８において、計算途中の電場および磁場の空間分布を共通部５６の電磁場として記憶部２０（ファイル）に保存する。

次に、ステップＳ１１０において、計算条件として、先に、ステップＳ１００で定義された計算モデルｉ（＝１，２，３）の屈折率分布を設定し、ステップＳ１０８で保存された共通部５６の電磁場の空間分布を読出／書込制御部２２によって記憶部２０（ファイル）から読み込む。

続いて、ステップＳ１１２において、ＦＤＴＤ法計算ステップ数を、ステップＳ１０４で設定した共通部５６の伝搬ステップ数Ｎ（＝１５３ｓｔｅｐ）に設定し、電磁場計算部１８によって、Ｎ目以降からの解析領域５０の電磁場計算を行ない、解析領域５０内の伝播する電磁場の時間的、空間的変化を計算していく。
なお、図示されてはいないが、計算モデル１の場合には、ステップＳ１１０をスルーして、ステップＳ１１２に移り、ステップＳ１０６における共通部５６の伝搬ステップ数Ｎ（＝１５３ｓｔｅｐ）の計算段階での計算の中断後、ステップＳ１１２における共通部５６の電磁場の記憶部２０（ファイル）への保存完了後に、もしくは保存中に、中断されていた電磁場計算部１８による電磁場の計算を、ステップＳ１１２において、Ｎ（＝１５３ｓｔｅｐ）以降から再開しても良い。

続いて、ステップＳ１１４において、計算後処理として、こうして計算された電磁場の時間、空間分布に基づいて、解析目的に応じた物理量を計算する。また、計算結果出力部２４において、計算結果を、表示デバイスのディスプレー上にグラフィック表示したり、プリンタにハードコピー出力する。
こうして、計算モデルｉについての電磁場計算が終了する。

次に、ステップＳ１１６において、全計算モデルの計算が終了したか否かを判断する。図示例では、ｉ＝３（または３以上）であるか否かを判断する。全計算モデルの計算が終了していれば（図示例では、ｉが３以上であれば）、Ｙｅｓとして、ステップＳ１１８（ＥＮＤ）に抜け出し、電磁場解析を終了する。一方、全計算モデルの計算が終了していなければ（図示例では、ｉが３未満であれば）、Ｎｏ：ｉ＝ｉ＋１として、ステップＳ１１０に戻る。

続いて、ステップＳ１１０において、計算条件として、先に、ステップＳ１００で定義された計算モデルｉ＋１（図示例では、２または３）の屈折率分布を設定し、記憶部２０に保存されている共通部５６の電磁場の空間分布を読出／書込制御部２２によって記憶部２０（ファイル）から読み込む。
以下、同様にして、ステップＳ１１６において全計算モデルの計算が終了したと判断されるまで、上記ステップＳ１１０〜Ｓ１１６を繰り返し、終了していればステップＳ１１８に抜け出して、電磁場解析を終了する。

このようにして、本実施例において、複数の計算モデルｉ（＝１，２，３）について行われた計算結果を、伝播ステップ数が５０（Ｎ＝５０ｓｔｅｐ）毎に、図８（ａ），（ｂ）および（ｃ）に示す。なお、図８（ｄ）には、計算モデルの解析領域のサイズを具体的に示す。
これらの図から、いずれの計算モデルにおいても、初期状態（０伝播ステップから１５０（正確には、１５３）伝播ステップ目までは、入射光が、解析領域５０の中心に位置する球形水滴粒子５２に到達しておらず、共通部５６を伝播しているので、共通部５６の中での電磁場は計算モデルによらず、同じである。

したがって、本実施例においては、共通部５６を伝播している初期状態（０伝播ステップから１５０（正確には、１５３）伝播ステップ目までの計算は、計算モデル１のみで行い、計算モデル２および３では省略できるので、図８（ｂ）および（ｃ）に示す初期状態（０伝播ステップから１５０（正確には、１５３）伝播ステップ目までの実線で囲まれた計算省略部分の電磁場の計算を省略できる。
その結果、図８に示す初期状態から終了（０伝播ステップから６００伝播ステップ目）まで計算を行う例においては、２つの計算モデルについて（０伝播ステップから１５３伝播ステップ目）までの計算を省略できるので、約１７％の計算量を削減できる。
計算量削減率＝（１５３×２）／（６００×３）×１００＝１７％

こうして、本発明においては、複数の計算モデルｉ（＝１，２，３）において、計算開始から、励震部付近の共通部を伝搬している間は、共通部の中での、電磁場は計算モデルによらず同じであるので、共通部の計算は、１つの計算モデル（ｉ＝１）のみにおいて行い、残りの計算モデルにおいては、共通部の計算を繰返さないことで、計算時間を短縮することができる。

上述したように、本発明は、同一の計算条件の下で、複数の異なる計算モデルに対する計算を連続して行なう場合の計算量の削減による計算の効率化についての技術であり、上述した光学、通信、半導体、検出素子、構造物などの技術分野に適用可能なことはもちろんであるが、例えば、計算モデル中のマイクロレンズなどにおいて、レンズの曲率、膜厚などを連続的に変化させて、最適条件を見出す場合などに行なわれる計算に、特に、有効である。
上述した例では、予め、屈折率分布が定義された複数の計算モデルについて、説明を行なったが、本発明はこれに限定されず、計算結果に応じて、屈折率分布を目的に合うように、自動的に修正していく、自動最適化の計算の場合にも適用できる。この場合には、屈折率分布が不変となる領域を、予め、先見情報として、設定しておくことで、同様の効率化を実現することができる。このような場合にも、類似した計算モデルに対する計算を、大量に反復して行なう必要があり、本発明が特に有効となる。

なお、上述した波動場解析方法は、予め用意されたプログラムをパーソナルコンピュータやワークステーションなどのコンピュータで実行することによっても実現されるものである。このプログラムは、波動場解析方法をコンピュータに実行させるプログラムであって、ハードディスク、フロッピーディスク、ＣＤ−ＲＯＭ、ＭＯ、ＤＶＤ等のコンピュータで読み取り可能な記録媒体に記録され、コンピュータによって記録媒体から読み出されることによって実行される。またこのプログラムは、上記記録媒体を介して、インターネット等のネットワークを介して配布することができる。

以上、本発明に係る波動場解析方法および装置ならびに波動場解析方法をコンピュータに実行させるプログラムを記録したコンピュータに読み取り可能な記録媒体について詳細に説明したが、本発明は、以上の実施形態に限定されるものではなく、本発明の要旨を逸脱しない範囲において、各種の改良や変更を行ってもよい。

本発明の波動場解析装置が適用される電磁場解析装置の機能的構成の一実施例を示すブロック図である。 本発明に適用されるＦＤＴＤ法を説明する説明図である。 本発明に適用されるＦＤＴＤ法よる磁場および電場の計算を説明する説明図である。 図１に示す電磁場解析装置のハードウエア構成の一実施例を示すブロック図である。 本発明の波動場解析方法が適用される電磁場解析方法の一例を示すフローチャートである。 （ａ），（ｂ）および（ｃ）は、それぞれ、図５に示す電磁場解析方法が適用される、異なる屈折率分布を持つ解析領域からなる計算モデルを示す模式図である。 （ａ），（ｂ）および（ｃ）は、それぞれ、複数の計算モデルの共通部の抽出を説明するための説明図である。 （ａ），（ｂ）および（ｃ）は、それぞれ、図６（ａ），（ｂ）および（ｃ）に示す計算モデルについて行われた計算結果を、所定伝播ステップ数毎に示す図であり、（ｄ）は、それらの計算モデルの具体的な寸法を説明する説明図である。 従来の電磁場解析方法を示すフローチャートである。

符号の説明

１０，３０ 電磁場解析装置
１２ 屈折率分布定義部
１４ 計算条件入力部
１６ 共通部抽出部
１８ 電磁場計算部
２０ 記憶部
２２ 書込／読出（制御）部
２４ 計算結果出力部
３２ 入力装置
３４ 読取／書込装置
３６ 演算処理装置（ＣＰＵ）
３８ メモリ
４０ 表示装置
４２ 記録装置（プリンタ）
４４ バス
５０，５０ａ，５０ｂ，５０ｃ 解析領域
５２，５２ａ，５２ｂ，５２ｃ 球形（水滴）粒子
５４ 励震部（波源）
５６ 共通部

Claims (13)

1. 解析領域内を波動が伝播する伝播状態を解析する波動場解析方法であって、
   同一の前記解析領域内の前記波動の伝播状態が異なる複数の計算モデルにおいて、初期状態から前記波動の伝播状態が共通である共通部を抽出し、
   抽出された前記共通部について、前記波動が伝播する前記解析領域内の波動場の時間変化の計算を初期状態から開始し、
   計算開始後、前記共通部において、前記解析領域内の波動場の時間変化を連続して計算し、
   前記共通部内の所定の途中の計算段階において、前記計算を停止して、計算された前記波動場の状態を記憶手段に保存し、
   前記複数の計算モデルの各モデルに対する、前記解析領域内の前記波動場の時間変化の計算を開始するに際し、前記記憶手段に保存された前記共通部内の前記所定の計算段階の前記波動場の状態を読み出し、
   前記記憶手段から読み出された前記共通部の前記所定の途中の計算段階の前記波動場の状態から、各モデルに対する、前記解析領域内の前記波動場の時間変化の計算を連続して行うことを特徴とする波動場解析方法。
2. 前記共通部内の所定の途中の計算段階の前記波動場の状態は、前記共通部の最後の計算段階の前記波動場の状態である請求項１に記載の波動場解析方法。
3. 前記波動は、前記解析領域内に設定された励震部において発生され、前記共通部は、前記励震部近傍の前記波動を伝播する媒質が共通である領域である請求項１または２に記載の波動場解析方法。
4. 前記波動は、電磁波であり、前記波動場は、電磁場である請求項１〜３のいずれかに記載の波動場解析方法。
5. 前記電磁波は、前記解析領域内に設定された励震部において発生され、
   前記共通部の抽出は、前記解析領域内に前記初期状態として初期の電磁場が存在せず、前記励震部近傍の屈折率分布が共通である領域を抽出することであり、
   前記共通部内の所定の途中の計算段階は、電磁波が前記共通部内を伝搬する時間に設定される請求項４に記載の波動場解析方法。
6. 前記複数の計算モデルが、前記電磁場の計算結果に基いて修正される請求項４または５に記載の波動場解析方法。
7. 解析領域内に波動が伝播する伝播状態を解析する波動場解析装置であって、
   前記波動が伝播する前記解析領域内の波動場の時間変化を初期状態から連続して計算する波動場計算手段と、
   同一の前記解析領域内の前記波動の伝播状態が異なる複数の計算モデルにおいて、前記初期状態から前記波動の伝播状態が共通である共通部を抽出する共通部抽出手段と、
   前記共通部抽出手段によって抽出された前記共通部について、前記波動場計算手段によって前記波動場の時間変化を前記初期状態から連続して計算することにより得られた、前記共通部内の所定の途中の計算段階の前記波動場の状態を保存する記憶手段と、
   前記記憶手段に保存された前記共通部内の前記所定の途中の計算段階の前記波動場の状態を読み出す読出手段とを有し、
   前記読出手段によって前記記憶手段から読み出された前記共通部の前記所定の途中の計算段階の前記波動場の状態から、前記波動場計算手段によって、前記複数の計算モデルの各モデルに対する、前記解析領域内の前記波動場の時間変化の計算を開始することを特徴とする波動場解析装置。
8. 前記共通部内の所定の途中の計算段階の前記波動場の状態は、前記共通部の最後の計算段階の前記波動場の状態である請求項７に記載の波動場解析装置。
9. さらに、前記解析領域内に設定され、前記波動を発生させる励震部を有し、
   前記共通部抽出手段は、前記共通部として、前記励震部近傍の前記波動を伝播する媒質が共通である領域を抽出する請求項７または８に記載の波動場解析装置。
10. 前記波動は、電磁波であり、前記波動場は、電磁場である請求項７〜９のいずれかに記載の波動場解析装置。
11. さらに、前記解析領域内に設定され、前記電磁波を発生させる励震部を有し、
    前記共通部抽出手段は、前記解析領域内に前記初期状態として、初期の電磁場が存在せず、前記励震部近傍の屈折率分布が共通である領域を抽出するものであり、
    前記共通部内の所定の途中の計算段階は、電磁波が前記共通部内を伝搬する時間に設定される請求項１０に記載の波動場解析装置。
12. 前記複数の計算モデルが、前記電磁場の計算結果に基いて修正される請求項１０または１１に記載の波動場解析装置。
13. 請求項１〜６のいずれかに記載の波動場解析方法をコンピュータに実行させるプログラムを記録したことを特徴とするコンピュータに読み取り可能な記録媒体。

Images