(57)【要約】【課題】本発明は、キャリア同士の相互作用のシミュレーションの処理負荷を軽減するキャリアシミュレーション装置を提供することを課題とする。【解決手段】キャリア同士の衝突による相互作用をシミュレーションするキャリアシミュレーション装置であって、エネルギと時間を変数とするキャリアの分布関数を、キャリア密度及びキャリア１個当たりの衝突時間間隔に基づく関数で定義し、変数であるエネルギと時間を離散値とし、キャリアの分布関数の値を有限個の数値の集合としてキャリア相互作用をシミュレーションすることを特徴とし、任意のエネルギをそれぞれ有するキャリア同士の衝突による相互作用によって新たなエネルギをそれぞれ有する各キャリアが生じる確率の左右非対称な分布を左右対称と近似するとともにキャリアが禁制帯へ遷移することを排除する関数を用いてキャリアの分布関数を定義すると好適である。【選択図】図８
2010177261.tif 000002
【特許請求の範囲】
【請求項１】
キャリア同士の衝突による相互作用をシミュレーションするキャリアシミュレーション装置であって、
エネルギと時間を変数とするキャリアの分布関数を、キャリア密度及びキャリア１個当たりの衝突時間間隔に基づく関数で定義し、
変数であるエネルギと時間を離散値とし、キャリアの分布関数の値を有限個の数値の集合としてキャリア相互作用をシミュレーションすることを特徴とするキャリアシミュレーション装置。
【請求項２】
任意のエネルギをそれぞれ有するキャリア同士の衝突による相互作用によって新たなエネルギをそれぞれ有する各キャリアが生じる確率の左右非対称な分布を左右対称と近似してキャリアの分布関数を定義することを特徴とする請求項１に記載するキャリアシミュレーション装置。
【請求項３】
キャリアが禁制帯へ遷移することを排除する関数を用いてキャリアの分布関数を定義することを特徴とする請求項２に記載するキャリアシミュレーション装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、キャリアシミュレーション装置に関する。
【背景技術】
【０００２】
太陽電池などのデバイスを設計する上でそのデバイスの動作を予測するためには、デバイスにおけるキャリアをシミュレーションすることが非常に重要となる。特許文献１には、半導体素子内部のキャリアの濃度分布、エネルギ分布、速度分布などの予測において、多重積分を回避するために、乱数を使って個々のキャリアの動きをシミュレーションするモンテカルロ法を使用することが記載されている。
【先行技術文献】
【特許文献】
【０００３】
【特許文献１】特開平８−１７８８５号公報
【特許文献２】特開２０００−１１０２７号公報
【発明の概要】
【発明が解決しようとする課題】
【０００４】
従来のキャリアについてのシミュレーションでは、新たに発生したキャリアが速やかにエネルギを失うことを前提とした物理モデルに基づいてシミュレーションを行っている。しかし、太陽電池などでは高効率化を図るために、エネルギを失わないで動作するキャリア（ホットキャリア）を利用することが検討されている。このような高いエネルギ状態を保持するキャリアの場合、キャリア同士が衝突して相互作用し、各キャリアのエネルギ状態が変化する。このようなキャリアの相互作用についての実用的なシミュレーション技術は、開発されていない。また、上記したようなモンテカルロ法を使用した場合でも、多数のキャリアについてのシミュレーションを行わないと平均的なキャリアの動作が判らないので、処理負荷やメモリが非常に膨大となる。
【０００５】
そこで、本発明は、キャリア同士の相互作用のシミュレーションの処理負荷を軽減するキャリアシミュレーション装置を提供することを課題とする。
【課題を解決するための手段】
【０００６】
本発明に係るキャリアシミュレーション装置は、キャリア同士の衝突による相互作用をシミュレーションするキャリアシミュレーション装置であって、エネルギと時間を変数とするキャリアの分布関数を、キャリア密度及びキャリア１個当たりの衝突時間間隔に基づく関数で定義し、変数であるエネルギと時間を離散値とし、キャリアの分布関数の値を有限個の数値の集合としてキャリア相互作用をシミュレーションすることを特徴とする。
【０００７】
このキャリアシミュレーション装置では、キャリアとキャリアとの衝突による相互作用によって変化するキャリアの分布をシミュレーションするために、キャリア密度とキャリア１個当たりの衝突時間間隔を用いてキャリアの分布関数を定義する。そして、キャリアシミュレーション装置では、このキャリアの分布関数にエネルギの離散値と時間の離散値を順次代入し、キャリアの分布関数の値を有限個の数値の集合として出力する。これによって、キャリアの分布の変化を、各時間（離散値）における各エネルギ（離散値）でのキャリアの分布として得ることができる。このように、キャリアシミュレーション装置では、キャリアの分布関数をエネルギと時間の離散値による有限個の数値集合として処理することにより、処理負荷を軽減することができる。また、キャリアシミュレーション装置では、キャリアの密度と衝突時間間隔を考慮したキャリアの分布関数と定義することにより、キャリア同士の相互作用についてのシミュレーション精度を確保することできる。
【０００８】
本発明の上記キャリアシミュレーション装置では、任意のエネルギをそれぞれ有するキャリア同士の衝突による相互作用によって新たなエネルギをそれぞれ有する各キャリアが生じる確率の左右非対称な分布を左右対称と近似してキャリアの分布関数を定義すると好適である。
【０００９】
任意のエネルギをそれぞれ有するキャリア同士が衝突して新たなエネルギをそれぞれ有する各キャリアが生じる過程における衝突前後の各エネルギの全ての組み合わせについてのキャリアが生じる確率を考慮してキャリアの分布状態のシミュレーションを行う場合、組み合わせが膨大となり、処理負荷が非常に大きくなる。そこで、このキャリアシミュレーション装置では、キャリア同士の衝突による相互作用によって新たなエネルギを有する各キャリアが生じる確率の左右非対称な分布を左右対称と近似してキャリアの分布関数を定義することにより、キャリアの分布関数を簡単化でき（例えば、キャリアの分布関数を多重積分計算式から１回程度の積分計算式に変換できる）、処理負荷を軽減することができる。
【００１０】
本発明の上記キャリアシミュレーション装置では、キャリアが禁制帯へ遷移することを排除する関数を用いてキャリアの分布関数を定義すると好適である。
【００１１】
キャリア同士の衝突による相互作用によって新たなエネルギを有する各キャリアが生じる確率の左右非対称な分布を左右対称と近似した場合、実際の物理現象（例えば、禁制帯へキャリアの分布が遷移しないこと）をシミュレーションできない虞がある。そこで、このキャリアシミュレーション装置では、キャリアが禁制帯へ遷移することを排除する関数を用いてキャリアの分布関数を定義することにより、キャリア同士の衝突による相互作用によって新たなエネルギを有する各キャリアが生じる確率の左右非対称な分布を左右対称と近似した場合でも実際の物理現象をシミュレーションすることができる。
【発明の効果】
【００１２】
本発明は、キャリアの分布関数をエネルギと時間の離散値による有限個の数値集合として処理することにより、処理負荷を軽減することができる。
【図面の簡単な説明】
【００１３】
【図１】本実施の形態に係るホットキャリアシミュレーション装置の構成図である。
【図２】一般的なキャリアの分布の時間変化の一例を示す図である。
【図３】エネルギｘ−ｙ平面上でのΩ関数を示す図である。
【図４】エネルギｘ軸でのΩ’関数を示す図である。
【図５】左右非対称なΩ’’関数を左右対称と近似した場合のキャリアの分布の時間変化の一例である。
【図６】Ｎ（Ｅ）関数の一例を示す図である。
【図７】Ｎ（Ｅ）関数によるキャリア分布の変化の一例を示す図である。
【図８】ホットキャリアシミュレーション処理の流れを示すフローチャートである。
【発明を実施するための形態】
【００１４】
以下、図面を参照して、本発明に係るキャリアシミュレーション装置の実施の形態を説明する。
【００１５】
本実施の形態では、本発明に係るキャリアシミュレーション装置を、パーソナルコンピュータ（以下、パソコンと記載）上に構成されるホットキャリアシミュレーション装置に適用する。本実施の形態に係るホットキャリアシミュレーション装置は、高効率なホットキャリア太陽電池の動作のシミュレーション（特に、キャリアのシミュレーション）を行う。
【００１６】
図１〜図７を参照して、本実施の形態に係るホットキャリアシミュレーション装置１について説明する。図１は、本実施の形態に係るホットキャリアシミュレーション装置の構成図である。図２は、一般的なキャリアの分布の時間変化の一例を示す図である。図３は、エネルギｘ−ｙ平面上でのΩ関数を示す図である。図４は、エネルギｘ軸でのΩ’関数を示す図である。図５は、左右非対称なΩ’’関数を左右対称と近似した場合のキャリアの分布の時間変化の一例である。図６は、Ｎ（Ｅ）関数の一例を示す図である。図７は、Ｎ（Ｅ）関数によるキャリア分布の変化の一例を示す図である。
【００１７】
ホットキャリアシミュレーション装置１は、記憶装置（ハードディスクなど）に格納されているホットキャリアシミュレーション用のプログラムをパソコン本体１０で実行することによって構成され、そのシミュレーション結果をディスプレイ１１に表示したりあるいは記憶装置に格納する。この際、ユーザは、このシミュレーションに必要な情報をキーボード１２などを利用して入力し、シミュレーション結果をディスプレイ１１などで確認する。
【００１８】
ホットキャリアシミュレーション装置１における処理について具体的に説明する前に、ホットキャリアシミュレーション方法について説明する。
【００１９】
一般的なキャリアは、新たに発生すると、速やかにエネルギを失ってゆき、高いエネルギ状態を保持できない。図２には、各グラフにおいて横軸をエネルギＥ、縦軸をキャリアの分布関数ｆ（Ｅ）の値とし、キャリアのエネルギ分布の時間変化を示している。図２から判るように、一般的なキャリアの場合、発生から時間が経過すると、エネルギが小さくなり、キャリアのエネルギ分布が禁制帯にへばりつくようになる。このようなキャリアを用いて太陽電池を構成した場合、キャリアのエネルギが小さく、移動が遅いので、発電効率が低い。そこで、発電効率を向上させるために、エネルギを失わず、高いエネルギ状態を保持することができるホットキャリアを用いて太陽電池を構成する。
【００２０】
ホットキャリアを用いた太陽電池の動作としては、以下の３つの特徴がある。１つ目の特徴は、キャリアがエネルギを失わずに、高いエネルギから低いエネルギにいたる広いエネルギ範囲にわたってキャリアが分布していることである。２つ目の特徴は、特定の狭い範囲のエネルギを持つキャリアを抜き取ることである。３つ目の特徴は、キャリアとキャリアとが衝突によって相互作用して各キャリアのエネルギが変化し、キャリアが抜き取られたことによるキャリア分布の穴を埋めることである。
【００２１】
このような太陽電池におけるキャリアについてのシミュレーションを行うには、キャリアの分布関数ｆ（Ｅ：エネルギ，ｔ：時間）を定義し、物理学におけるボルツマン方程式に基づいてそのキャリアの変化を数値計算すればよい。式（１）は、キャリアの分布関数ｆ（Ｅ，ｔ）を用いたボルツマン方程式である。
【数１】
2010177261.tif 000003
【００２２】
式（１）の右辺の第１項は、キャリアの消滅に相当する項であり、Ｔｌ（Ｅ）がキャリアが消滅する平均時間を表わす時定数である。Ｔｌは、エネルギＥに応じて変化する時定数としているが、固定の時定数でもよい。式（１）の右辺の第２項は、キャリアの減衰に相当する項であり、α（Ｅ）が単位時間Ｔｐ（Ｅ）当たりに減少する平均エネルギＥｐ（Ｅ）を表わす。式（１）の右辺の第３項は、キャリア移動層からエネルギとして取り出されるキャリアの消滅に相当する項であり、Ｔｅ（Ｅ）がキャリアが取り出される平均時間を表わす時定数である。式（１）の右辺の第４項は、光を受けて新たに生まれるキャリアの供給に相当する項であり、Ｔｓ（Ｅ）がキャリアが新たに生まれる平均時間を表わす時定数である。Ｔｓは、エネルギＥに応じて変化する時定数としているが、固定の時定数でもよい。式（１）の右辺の第５項は、キャリア同士の衝突による相互作用に相当する項である。なお、式（１）の右辺の第１項〜第４項は、従来の手法によって簡単に計算でき、処理負荷も軽い。
【００２３】
以下で、式（１）の右辺の第５項のキャリア同士の衝突項について考える。キャリア同士が衝突した場合、ａのエネルギを持ったキャリアとｂのエネルギを持ったキャリアが相互作用し、ｘのエネルギを持ったキャリアとｙのエネルギを持ったキャリアが新たに生じる。キャリア同士の衝突項では、この衝突過程における４つのエネルギａ，ｂ，ｘ，ｙの全ての組み合わせについてキャリアが生じる確率を考慮しながらキャリアのエネルギ分布の変化を考える必要がある。このキャリア同士の相互作用によるエネルギ変化の分布を示す確率関数をΩ（ａ，ｂ，ｘ，ｙ）とすると、式（１）の右辺の第５項のキャリア同士の衝突項は式（２）となる。式（２）におけるＴｃは、キャリア１個当たりの衝突の平均時間間隔である。
【数２】
2010177261.tif 000004
【００２４】
式（２）から判るように、４つのエネルギａ，ｂ，ｘ，ｙの全ての組み合わせについての確率を考慮する場合、組み合わせが膨大であり、４階積分計算式となる。式（２）をそのまま計算した場合、処理負荷が非常に大きくなり、記憶容量も大容量必要となる。したがって、実用的なシミュレーションを行うためには、この式（２）を、ある程度のシミュレーション精度を確保した上で、簡略化し、処理負荷や記録容量を軽減する必要がある。
【００２５】
図３には、横軸を一方のキャリアのエネルギとし、縦軸を他方のキャリアのエネルギとし、エネルギｘ−ｙ平面上での確率関数Ω（ａ，ｂ，ｘ，ｙ）を示している。エネルギａを持つキャリアとエネルギｂを持つキャリアとの衝突はその逆のエネルギを持つ場合も同じであり、その２つのキャリアを区別することはできない。したがって、図３に示すように、確率関数Ωと線対称で同じ確率関数Ωが存在する。
【００２６】
式（２）におけるＴｃがキャリア１個当たりの衝突平均時間間隔であるので、キャリアの分布の確率関数Ω（ａ，ｂ，ｘ，ｙ）を積分すると１のオーダーとなる。ところが、エネルギの保存則からａ＋ｂ＝ｘ＋ｙであるので、式（３）が成立する。式（３）における関数δ（ｘ）は、０で無限大になる関数であり、式（３’）を満たす関数である。
【数３】
2010177261.tif 000005
【００２７】
上記したようにエネルギａを持つキャリアとエネルギｂを持つキャリアとの衝突ではその２つのキャリアを区別できないので、式（３）のΩ’（ａ，ｂ，ｘ）はｘに対して対称な２つの山が重なってできている。ここで、物理的な近似を導入する。確率関数Ω’（ａ，ｂ，ｘ）は、キャリアのエネルギレベルに依存しないというものである。図４は、横軸をエネルギｘとし、縦軸をΩ’（ａ，ｂ，ｘ）とし、エネルギｘ軸上のΩ’を示している。実際には、例えば、低エネルギ（移動が遅い）のキャリア同士よりも高エネルギ（移動が速い）のキャリア同士の衝突の方の確率が高い、あるいは、高エネルギのキャリアと低エネルギのキャリアの相互作用では互いにエネルギを平均化するような反応の確率が高い。したがって、このような確率の分布は、左右対称とはならず、図４の２つの山（Ω’’（ｘ−ａ）、Ω’（ｂ−ｘ））でそれぞれ示すように左右非対称になっている。
【００２８】
ここで、式（２）を簡単化するために、この確率の左右非対称な分布を無視し、左右対称と近似する。この近似は、キャリアの伝導帯のエネルギ範囲を−∞〜∞であるとみなす結果にもなる。したがって、この近似は、物理的にはかなり大きな制約であるが、数学的にはドラスティックな扱いを容易にすることができる。この近似により、Ω’（ａ，ｂ，ｘ）を式（４）と置くことができる。
【数４】
2010177261.tif 000006
【００２９】
式（２）に式（３）と式（４）を代入すると、式（５）のように変形できる。
【数５】
2010177261.tif 000007
【００３０】
式（５）の右辺の積分における第１項、第３項及び第４項は、ｙで積分すると、単にδ関数が消えるだけになる。また、式（５）の右辺の積分における第２項は、ｙを含んでいないので、ｙで積分すると何も変化が無く、さらに、ｘで積分するとｘ＝ａ＋ｂ−Ｅを各項に代入した形になる。以上の演算を行うと、式（５）を式（６）に変形できる。
【数６】
2010177261.tif 000008
【００３１】
式（６）の右辺の積分における第１項と第２項は、同じである。また、式（６）の右辺の積分における第３項では、ｂが単に積分変数だから、ｂをａと書き換えてもよい。したがって、式（６）を式（７）に変形できる。
【数７】
2010177261.tif 000009
【００３２】
式（７）を並び替え、式（８）に変形する。但し、式（７）の右辺の積分における第２項のａとｂの表記は入れ替えてある。
【数８】
2010177261.tif 000010
【００３３】
式（８）を更に並び替え、式（９）に変形する。
【数９】
2010177261.tif 000011
【００３４】
式（９）の右辺の第２項の積分は、−∞〜∞の積分なので、ｘ’＝ｘ−ａと置いてもよい。同様に、式（９）の右辺の第３項の積分も、−∞〜∞の積分なので、ｘ’’＝ｘ−Ｅと置いてもよい。したがって、式（９）を式（１０）に変形できる。
【数１０】
2010177261.tif 000012
【００３５】
式（１０）の右辺の第２項のｘ’による積分では、ｆ（ａ）を外に出せる。同様に、式（１０）の右辺の第３項のｘ’’による積分でも、ｆ（ａ）を外に出せる。したがって、式（１０）を式（１１）に変形できる。
【数１１】
2010177261.tif 000013
【００３６】
この式（１１）を整理すると、式（１２）に変形できる。
【数１２】
2010177261.tif 000014
【００３７】
ここで、式（１３）、式（１４）、式（１５）で示す定義を行う。
【数１３】
2010177261.tif 000015
【００３８】
式（１２）に、式（１３）、式（１４）、式（１５）の定義を適用すると、式（１６）に変形できる。但し、式（１６）ではａからｘへ記号変更している。
【数１４】
2010177261.tif 000016
【００３９】
式（１６）は、式（４）の近似を除いて最も一般化したキャリア同士の相互作用の表現であり、式（４）の近似に基づいて式（２）で示す衝突項を変形したものである。式（１６）におけるρは、式（１５）の定義から明らかなように、キャリアの密度である。式（１５）から判るように、ρは、ｆ（Ｅ）に依存する量であるが、計算では定数のように扱って差し支えない。したがって、式（１６）は、実質的には線形の１階積分であり、式（２）の非線形の４階積分に比べるとはるかに扱い易い。そのため、計算負荷が非常に少なくなり、計算に必要なメモリも少なくてすむ。
【００４０】
式（１６）のモデルを構築するにあたり、式（４）の近似を行った。実際には、式（４）の近似には、上記したように影響がある。特に、量子力学効果によって、キャリアの分布が禁止されたエネルギレベル（禁制帯）への遷移が無いことや他の粒子で占められた状態への遷移がパウリの排他律によって禁止されるという物理現象は重要であるが、このような物理現象を式（４）の近似によって織り込めない。式（１６）のままでシミュレーションを行うと、初期条件として禁制帯を設けたキャリア分布を与えても、図５に示すように、キャリアが禁制帯に遷移し、禁制帯へのキャリア分布がどんどん増えるという結果になってしまい、実際の物理現象を再現できない。そこで、禁制帯へのキャリア分布の遷移を防止するために、関数Ｎ（Ｅ）を導入し、式（１４）の定義を式（１７）に変更する。そして、式（１６）に、式（１７）の定義を適用し、式（１８）に変形する。
【数１５】
2010177261.tif 000017
【００４１】
式（１７）におけるＮ（Ｅ）は、状態密度を考慮した関数であり、エネルギＥの準位への遷移し易さを表わす関数である。Ｎ（Ｅ）は、禁制帯のエネルギ準位ではゼロとする関数であり、これによって禁制帯へのキャリア分布の遷移を無くすことをモデル化する。Ｎ（Ｅ）は、例えば、図６に示すような関数である。また、Ｎ（Ｅ）としては、最も簡単な関数として、禁制帯では０とし、伝導帯では１とする関数である。このように、Ｎ（Ｅ）は、禁制帯で実質的にゼロである関数である。
【００４２】
式（１７）と等価な定義として、式（１９）の定義としてもよい。この式（１９）では、ｘ＞ＥｃでＮ＝１であり、ｘ＜ＥｃでＮ＝０としたものと同じである。また、式（１７）と等価な定義として、式（２０）の定義としてもよい。この式（２０）では、Ｎ（Ｅ）＝１／Ｍ（Ｅ）と表現を変えただけである。
【数１６】
2010177261.tif 000018
【００４３】
なお、上記したパウリの排他律の影響は、モデル化できないが、伝導帯のキャリア密度は元々パウリの排他律が影響するような密度ではないため、問題にならない。
【００４４】
この式（１７）、式（１８）のモデルにおける粒子数保存則について考える。粒子（キャリア）の全数は、キャリアの消滅（式（１）の右辺の第１項）、キャリア移動層からエネルギとして取り出されるキャリアの消滅（式（１）の右辺の第３項）及び光を受けて新たに生まれるキャリアの供給（式（１）の右辺の第４項）のため、保存しているわけではない。しかし、キャリア同士の衝突による相互作用によって、粒子数が増減してはならない。したがって、式（２１）の条件を満たしていることが必要である。
【数１７】
2010177261.tif 000019
【００４５】
ここで、式（１７）、式（１８）が式（２１）を満たしていることを検証する。式（１７）、式（１８）により、式（２２）となる。
【数１８】
2010177261.tif 000020
【００４６】
ここで、式（２２）の右辺の積分におけるｘとＥの記号を入れ替えてもよく、また、Ω’’’（ｘ）は対称な関数であり、Ω’’’（ｘ）＝Ω’’’（−ｘ）なので、式（２２）の右辺の［］の中は、等式の引き算となり、０になる。したがって、式（１７）、式（１８）は、粒子数保存則を満たしている。
【００４７】
この式（１７）、式（１８）のモデルにおけるエネルギ保存則について考える。粒子数保存則と同様に、全エネルギについてもキャリア同士の衝突による相互作用によって増減してはならない。したがって、式（２３）の条件を満たしていることが必要である。
【数１９】
2010177261.tif 000021
【００４８】
ここで、式（１７）、式（１８）が式（２３）を満たしていることを検証する。式（１７）、式（１８）及びΩ’’’（ｘ）＝Ω’’’（−ｘ）、ｘとＥの記号を入れ替えにより、式（２４）となる。
【数２０】
2010177261.tif 000022
【００４９】
この式（２４）の右辺は、０にならない。したがって、式（１７）、式（１８）は、一般に、エネルギ保存則を満たしていない。エネルギ保存が壊れるのは、Ｎ（Ｅ）が一定でない領域である。図７には、横軸をエネルギＥとし、Ｎ（Ｅ）によるキャリア分布を示している。Ｎ（Ｅ）の変化している部位において、ｆ（Ｅ）が直線的変化であったとすると、Ω’’’（ｘ）が狭い関数であった場合、Ｎ（Ｅ）のグラフが凸のところではｆ（Ｅ）は増加し、Ｎ（Ｅ）のグラフが凹のところではｆ（Ｅ）は減少する。式（１８）でｆ（Ｅ）を一定とすると、式（１８）を式（２５）に変形できる。
【数２１】
2010177261.tif 000023
【００５０】
この式（２５）においてΩ’’’をほぼ定数とする。この場合、Ｎ（Ｅ）のグラフが凸のところでは式（２５）の値は正になるので、ｆ（Ｅ）は時間とともに増加する。一方、Ｎ（Ｅ）のグラフが凹のところでは式（２５）の値は負になるので、ｆ（Ｅ）は時間とともに減少する。これにより、ｆ（Ｅ）とΩ’’’が実際には変数である場合においても、上記の傾向があることを理解できる。
【００５１】
図７の場合、相対的に高エネルギのキャリアの分布が増え、低エネルギのキャリアの分布が減る分だけ、総エネルギが増加する。この式（１７）、式（１８）によるモデルでは、Ｎ（Ｅ）の値の小さいところからキャリアを追い出してＮ（Ｅ）とｆ（Ｅ）が同じ形のグラフにしようとし、グラフの形が同じになった後は変化が無く、エネルギ変化も起こらなくなる。
【００５２】
実際の解析では、Ｎ（Ｅ）が変化する領域は、伝導帯の端部に限られ、定常解を求めることが主なので、エネルギが保存しないという問題は起こり難く、影響は少ない。特に、Ｎ（Ｅ）が０か１のどちらかの値しか取らないような場合、Ｎ（Ｅ）＝０である領域ではｆ（Ｅ）も０なので、式（２４）の積分に影響せず、積分結果は０になって完全にエネルギ保存則は満たされる。
【００５３】
したがって、式（１７）、式（１８）によるモデルは、キャリア分布の禁制帯への遷移を排除するとともに、粒子数保存則を満たしており、エネルギ保存則への影響も少ないので、実質的なシミュレーション精度を十分に確保できる。また、式（１７）、式（１８）によるモデルは、線形の１階積分であるので、計算が簡便であり、処理負荷を軽減できる。
【００５４】
それでは、ホットキャリアシミュレーション装置１におけるホットキャリアのシミュレーション処理について説明する。ホットキャリアシミュレーション用のプログラムをパソコン本体１０で実行することにより、ホットキャリアシミュレーションのルーチンが実行される。なお、このホットキャリアシミュレーションのルーチンを、太陽電池の動作のシミュレーションを行うためのメインルーチンから呼び出されて実行されるサブルーチンとしてもよい。
【００５５】
関数ｆ（Ｅ，ｔ）は、原理的には無限の自由度を持つため、そのままパソコンなどで計算を行うことは不可能である。そこで、ホットキャリアシミュレーション装置１では、エネルギＥと時間ｔを適当な値に分割して離散値として取り扱い、関数ｆ（Ｅ、ｔ）を有限個の数値の集合として計算する。
【００５６】
実際のキャリアの現象では、エネルギＥは基底の０から４ｅＶ程度まで考慮すれば十分であり、その間を例えば２００分割程度すれば十分な精度でシミュレーションが可能である。したがって、ホットキャリアシミュレーション装置１では、エネルギＥの離散値として、Ｅｍａｘ＝４ｅＶ、Ｅｍｉｎ＝０、ｄＥ＝０．０２ｅＶとする。なお、Ｅｍａｘ、Ｅｍｉｎ、ｄＥについは、一例であり、これ以外の数値でもよい。
【００５７】
また、時間ｔに関しては、キャリアの現象は１００ｆｓｅｃ程度の時間オーダで現象が進行することが判っており、それより十分に小さい時間、例えば、ｄｆ＝５０ｆｓｅｃ程度とすれば十分な精度でシミュレーションが可能である。シミュレーションしたい時間スケールは、最大１０００ｐｓｅｃ程度であり、時間分割は２万分割である。したがって、ホットキャリアシミュレーション装置１では、時間ｔの離散値として、ｔｍａｘ＝１０００ｐｓｅｃ、ｔｍｉｎ＝０、ｄｔ＝５０ｆｓｅｃとする。なお、ｔｍａｘ、ｔｍｉｎ、ｄｔについは、一例であり、これ以外の数値でもよい。
【００５８】
なお、このような離散値を用いた場合でも多くの記憶容量が必要となるため、例えば、１０ｐｓｅｃ毎に計算結果を固定的な記憶装置に蓄えて記録し、その途中の値（Ｅに関して２００個の値）は揮発性の記憶装置に一時的に蓄え、最終的にｆ（Ｅ，ｔ）の値として２万個を固定的な記憶装置に格納する。
【００５９】
物理的なモデルとしてＴｃ、Ｎ（Ｅ）、Ω（Ｅ）、Ｒ（Ｅ）が必要であるが、Ｒ（Ｅ）はＮ（Ｅ）とΩ（Ｅ）から算出することができる。Ｎ（Ｅ）とΩ（Ｅ）は、ユーザによって予め定義（入力）され、ｆ（Ｅ，ｔ）と同様にＥに関して分割した値の集合として記憶装置に予め格納される。Ｔｃは、１個の値であり、ユーザによって予め定義（入力）され、記憶装置に予め格納される。ここでは、Ｒ’（Ｅ）をＲ（Ｅ）と改め、Ω’’’（Ｅ）をΩ（Ｅ）と改め、式（１７）を式（２６）、式（１８）を式（２７）として用いる。
【数２２】
2010177261.tif 000024
【００６０】
それでは、ホットキャリアシミュレーション処理を図８のフローチャートに沿って説明する。図８は、ホットキャリアシミュレーション処理の流れを示すフローチャートである。ここでは、ｘ，ｙは、積分のエネルギ変数である。ｄｆは、次の時間へのｆの変化値である。ｄｆｏは、ボルツマン方程式（式（１））におけるキャリアの衝突項以外のｆの変化値である。なお、フローチャート中で、ｆ（Ｅ，ｔ）は、Ｅとｔを分割した集合の中でＥとｔが一致する記憶値に対して、値の代入（記憶値の更新）及び値の呼び出し（記憶値の呼び出し）を指す。より実際的なプログラムでは、ｆ（Ｅ，ｔ）は整数の番号ｍ，ｎによって、例えば、Ｆ（ｍ、ｎ）という形で管理され、ｍ＝Ｅ／ｄＥ、ｎ＝ｔ／ｄｔとして取り扱う。
【００６１】
ルーチンが起動すると、ｆ（Ｅ，０）を初期化する（Ｓ１）。また、記憶装置に予め格納されているエネルギＥ毎のＮ（Ｅ）とΩ（Ｅ）を用いて、式（２６）により、エネルギＥ毎のＲ（Ｅ）を算出し、記憶装置に格納する（Ｓ２）。
【００６２】
初期化後、ループ毎に、ｄｆに０を代入し、ρに０を代入し、ｘにＥｍｉｎを代入し、ｙにＥｍｉｎを代入して初期化する（Ｓ３）。
【００６３】
次に、式（１５）に基づいて、ｆ（ｙ，ｔ）×ｄＥを算出し、その算出値をρの前回値に加算し、その加算値をρに代入する（Ｓ４）。また、ｙの今回値にｄＥを加算し、その加算値をｙに代入する（Ｓ５）。次に、ｙがＥｍａｘより小さいか否かを判定する（Ｓ６）。Ｓ６にてｙがＥｍａｘより小さいと判定した場合、Ｓ４の処理に戻り、次のエネルギの値についての処理を行う。
【００６４】
Ｓ６にてｙがＥｍａｘと等しくなったと判定した場合、式（１）に基づいて、既存のシミュレーションによってｄｆｏを算出する（Ｓ７）。
【００６５】
次に、計算済みのρ、記憶装置に予め格納されている今回のエネルギの値でのＮ（Ｅ）とΩ（Ｅ）及びＴｃを用いて、式（２７）に基づいて、ρ×ｆ（ｘ，ｔ）×Ω（ｘ−Ｅ）×Ｎ（ｘ）×ｄＥ／Ｔｃを算出し、その算出値をｄｆの前回値に加算し、その加算値をｄｆに代入する（Ｓ８）。また、ｘの今回値にｄＥを加算し、その加算値をｘに代入する（Ｓ９）。次に、ｘがＥｍａｘより小さいか否かを判定する（Ｓ１０）。Ｓ１０にてｘがＥｍａｘより小さいと判定した場合、Ｓ８の処理に戻り、次のエネルギの値についての処理を行う。
【００６６】
Ｓ１０にてｘがＥｍａｘと等しくなったと判定した場合、算出済みのＲ（Ｅ）とρを用いて、式（２７）に基づいて、ρ×Ｒ（Ｅ）×ｆ（Ｅ，ｔ）／Ｔｃを算出し、その算出値をＳ８〜Ｓ１０で算出したｄｆから減算し、その減算値を算出済みのｄｆｏに加算し、その加算値をｄｆに代入する（Ｓ１１）。
【００６７】
次に、Ｓ１１で算出したｄｆを用いて、ｄｆ×ｄｔを算出し、その算出値をｆ（Ｅ，ｔ）の前回値に加算し、その加算値をｆ（Ｅ，ｔ）に代入し、そのｆ（Ｅ，ｔ）を今回のエネルギＥと時間ｔのｆ（Ｅ，ｔ）の値として揮発性の記録装置に格納する（Ｓ１２）。
【００６８】
次に、Ｅの今回値にｄＥを加算し、その加算値をＥに代入する（Ｓ１４）。次に、ＥがＥｍａｘより小さいか否かを判定する（Ｓ１４）。Ｓ１４にてＥがＥｍａｘより小さいと判定した場合、Ｓ７の処理に戻り、次のエネルギの値についての処理を行う。
【００６９】
Ｓ１４にてＥがＥｍａｘと等しくなったと判定した場合、ｔの今回値にｄｔを加算し、その加算値をｔに代入する（Ｓ１５）。次に、ｔがｔｍａｘより小さいか否かを判定する（Ｓ１６）。Ｓ１６にてｔがｔｍａｘより小さいと判定した場合、Ｓ３の処理に戻り、次の時間の値についての処理を行う。
【００７０】
Ｓ１６にてｔがｔｍａｘと等しくなったと判定した場合、離散的なエネルギと時間の全ての組み合わせについてのｆ（Ｅ，ｔ）が算出され、この全ての組み合わせについてのｆ（Ｅ，ｔ）が記憶装置に格納されている。ホットキャリアシミュレーション装置１では、この記憶装置に格納されている時間毎の各エネルギに対するｆ（Ｅ，ｔ）の値を、数値表としてディスプレイ１１に表示したり、図２に示すようなキャリアのエネルギ分布グラフの時間変化としてディスプレイ１１に表示するなどして、ユーザに提供する。また、ホットキャリアシミュレーション装置１では、この記憶装置に格納されている時間毎の各エネルギに対するｆ（Ｅ，ｔ）の値を、太陽電池の動作のシミュレーションに出力する。
【００７１】
このホットキャリアシミュレーション装置１によれば、キャリアの衝突項として式（１７）、式（１８）による簡便なモデル式を構築し、式（１７）、式（１８）に基づいてエネルギＥと時間ｔの離散値で有限個のｆ（Ｅ，ｔ）の値を得ることにより、処理負荷を軽減でき、実用的なホットキャリアシミュレーションを行うことができる。さらに、ホットキャリアシミュレーション装置１では、式（１７）、式（１８）によるモデル式でシミュレーションを行うことにより、実質的なシミュレーション精度も十分に確保できる。
【００７２】
また、ホットキャリアシミュレーション装置１によれば、確率関数Ω’’の左右非対称な分布を左右対称と近似することにより、式（２）で示す非線形の４階積分式を式（１８）で示す線形の１階積分式まで変換でき、処理負荷を大幅に軽減することができる。また、ホットキャリアシミュレーション装置１によれば、関数Ｎ（Ｅ）を導入して最終的に式（１８）のモデル式を構築することにより、確率関数Ω’’の左右非対称な分布を左右対称と近似しても、キャリア分布の禁制帯への遷移を防止でき、シミュレーション精度を確保できる。
【００７３】
以上、本発明に係る実施の形態について説明したが、本発明は上記実施の形態に限定されることなく様々な形態で実施される。
【００７４】
例えば、本実施の形態ではホットキャリアシミュレーションのプログラムを実行することによってパソコン上に構成されるホットキャリアシミュレーション装置に適用したが、ホットキャリアシミュレーション装置の専用機として構成してもよいし、あるいは、ＣＤ−ＲＯＭなどの記憶媒体に格納されたプログラム、インタネットなどのネットワークを介して利用可能なプログラムなどに適用してもよい。
【００７５】
また、本実施の形態では太陽電池で利用するホットキャリアのシミュレーションに適用したが、太陽電池以外にも様々なものに適用可能（例えば、ビデオ撮像素子）であり、また、ホットキャリア以外にも適用可能（例えば、一般的なキャリア）であり、また、定常状態のキャリアのシミュレーションにも適用可能である。
【００７６】
また、本実施の形態ではキャリアの衝突項として式（１７）、式（１８）で示すモデル式を構築したが、このモデル式以外でもよい。
【００７７】
また、本実施の形態ではキャリアが禁制帯へ遷移することを排除するための関数としてＮ（Ｅ）を導入したが、禁制帯で実質的に０になるような関数であれば他のものでもよい。
【符号の説明】
【００７８】
１…ホットキャリアシミュレーション装置、１０…パーソナルコンピュータ本体、１１…ディスプレイ、１２…キーボード
【図１】
2010177261.tif 000025
【図２】
2010177261.tif 000026
【図３】
2010177261.tif 000027
【図４】
2010177261.tif 000028
【図５】
2010177261.tif 000029
【図６】
2010177261.tif 000030
【図７】
2010177261.tif 000031
【図８】
2010177261.tif 000032