2次元時系列差分法を用いたペロブスカイト太陽電池チュートリアル：時系列差分法とドリフト拡散法の連成シミュレーション

1. はじめに

このチュートリアルでは、テクスチャ処理が施されたペロブスカイト太陽電池に対する2次元有限差分時域法（FDTD）およびドリフト拡散法の連成シミュレーションを紹介する。FDTDソルバーによりデバイス内部の光電場および吸収特性が算出される。得られた空間的な光生成率のデータは電気的ソルバーへ直接渡され、そこでキャリア輸送、再結合、電荷抽出の過程が計算され、最終的に電流密度電圧特性（JV曲線）が求められる。

この太陽電池デバイスには??に示すような光閉じ込め用の表面構造が備わっている。この構造により前面表面での反射が抑制され、光が斜め方向へと屈折され、ペロブスカイト吸収層内での光路長が増大する。その結果、吸収層の厚みを増加させることなく光吸収効率が向上する。厚み増加は再結合を促進しキャリア収集率を低下させる可能性があるため、この構造は非常に有益である。

転送行列法のような一次元光学モデルは、平坦で横方向に均一な層構造に適している。しかし、テクスチャー加工された表面からの横方向散乱を完全に記述することはできない。本デバイスにおいては、電磁場がセル全体および厚み方向の両方で変化するため、空間メッシュ上でマクスウェル方程式を解くための二次元FDTDモデルが用いられる。

FDTD法により局所的な電界が得られる\(\mathbf{E}(x,y)\)。OghmaNanoはこれを基に吸収された光パワー密度を算出し、それを空間的な光生成率\(G(x,y)\)へと変換する\(\mathrm{m^{-3}\,s^{-1}}\)。その後、ドリフト拡散ソルバーが\(G(x,y)\)を用いて電圧を変化させながらキャリア密度、静電ポテンシャル、再結合現象、電流の流れ、および準フェルミ準位を計算する。

このような連成解析手法は極めて重要である。なぜなら、光吸収量が増加しても必ずしも光電流も同様に増加するわけではないからだ。再結合が活発であるか輸送特性が劣る領域で生成されたキャリアは効率的に収集されない可能性がある。したがって本シミュレーションでは、単なる吸収率だけでなく最終的な電気出力の観点から光閉じ込め構造の効果を評価しているのである。

このチュートリアルでは、完全な光電子シミュレーションを実行し、デバイスのメッシュ構造および光生成分布を確認するとともに、電荷キャリアの分布状態を調べ、光捕捉構造の有無における電流密度電圧特性を比較します。

2. 新規シミュレーションの作成

「新規シミュレーション」ウィンドウを開き、FDTD例のカテゴリを選択してください（??）。次にペロブスカイト用FDTD例を選択します（??）。シミュレーションを適切なフォルダに保存すると、??に示すメインインターフェースが表示されます。

3. メインウィンドウの見方

シミュレーションを読み込んだ後、??に示すようにメインウィンドウが表示されます。デバイス構造タブ内の3Dビューに構造物が描出されます。黒色の背景はシミュレーション空間、すなわちデバイスが存在する領域を表しています。

いくつかの重要な構成要素が確認できる。デバイス上部にある緑色の矢印群がFDTD光源を表している。光は上方から入射し、実際の太陽電池に日光が当たるのと同様にセル内部へと伝播する。光源の下には検出器面が存在し、光が通過する際の電磁場を測定する。この測定値は、後程デバイス内部の光生成分布を再構築するために不可欠となる。

構造体の縦方向に紫色の格子が存在する。これがFDTD用のメッシュであり、シミュレーションが実行される断面面を表している。つまり、デバイスの2次元断面を取り出して解析しているのである。セルの上部付近には三角形形状の光捕捉構造も見られるが、これは入射光を散乱させ吸収率を高めることを目的としている。その下には、透明電極、電子・正孔輸送層、ペロブスカイト吸収層、金属製裏面電極といった機能性ペロブスカイト層が積層されており、実際に動作する太陽電池を形成している。

なお、検出器は光そのものに干渉することはなく、単に電磁場を測定するだけである。

まもなくシミュレーションを実行し、光が構造体内部へ進入して光閉じ込め層で散乱された後、ペロブスカイト層で吸収される様子を観察します。この吸収作用によって電荷キャリアが生成され、電気モデルによって収集されるのです。

4. シミュレーションの実行

シミュレーションを保存し、構造の確認が済んだら「Run simulation (▶)」ボタンをクリックして計算を開始します。ウィンドウが「Terminal」タブに切り替わり、計算の進行状況を確認できるようになります。シミュレーションは複数の段階を経て実行され、その様子は??から??に示されています。また、計算結果は??に記載されています。

校正用シミュレーション。FDTDモジュールが起動すると（??）、まず最初に校正用シミュレーションが実行される。この段階においてOghmaNanoは一時的にデバイスを計算領域から除外し、光源パルスを空間内で伝播させる。これにより検出器において入射光場の基準値が記録される。つまり、デバイスが存在しない場合にどの程度の光が到達するかが得られるのだ。その後、プログラムは実際のシミュレーション（デバイスが存在する状態）の結果とこの基準値とを比較する。この比較処理によってOghmaNanoは反射光、透過光、そして実際に吸収された光を分離でき、結果としてセル内部の光生成分布を正確に算出することが可能となる。校正処理が実行中は、出力される各行の先頭にcalibrationという文字列が付加される。

本番の光学シミュレーション。キャリブレーションが完了すると（??）、出力内容からcalibrationという文字が消え、実際のデバイスに対するFDTDシミュレーションが開始され、電磁場が時間の経過とともに更新されていく。出力される進行状況のテキストが緑色で表示されている点に注目されたい。これは意図的な処置であり、本例では波長約550 nmの単一波長のみで解析が行われている。この波長は可視光の緑色領域に相当するため、OghmaNanoでは現在シミュレーション中の波長に合わせてタイムステップ情報の色を緑色に設定しているのである。

電気シミュレーションの実行。FDTD解析が完了すると、端末に経過時間が表示される。この例では光学的解析に約30秒を要した（??）。FDTD法は構造内での光の伝播を時刻ごとに明示的にシミュレートするため、本質的に計算速度が遅いという特徴がある。光学的解析が終了すると、OghmaNanoは自動的に電気シミュレーションへ移行し、FDTD法によって算出された光生成率分布を用いて電流密度電圧特性曲線を求める。端末には印加電圧が0 Vから順次増加していく様子とそれに伴う電流密度の値が表示される。電流密度が負の値を示す場合、これはセルが電流を生成していることを意味し、すなわち照明下において正常に機能する太陽電池として動作していることを表す。

計算が完了したら、出力タブ（??）に切り替えてシミュレーションで生成されたデータを確認します。次のセクションでは、最も重要な出力結果について説明します。

5. 結果の確認方法

計算が完了したら、??タブから出力結果を確認できます。ここでは特に情報量の多い3つの結果、すなわちデバイスメッシュ、光生成分布、および電流電圧特性に注目します。

デバイスメッシュ。まず device を開いて ?? に示される幾何形状を確認します。これにより、OghmaNanoが内部で構造をどのように表現しているかが明らかになります。すなわち、デバイス全体が多数の三角形に分割されているのです。この三角形メッシュこそがソルバーが実際に認識する構造であり、回転や拡大縮小を自由に行うことで、各層や光閉じ込め構造がどのように構築されているかを詳細に調べることが可能です。

光生成分布。次に G を開いて ?? に示されている光生成率を確認してください。これはFDTD法によって算出された、デバイス内部で電子・正孔対が生成される場所の空間分布図です。この値は \(\mathrm{m^{-3}\,s^{-1}}\) という単位、すなわち単位体積あたり1秒間に生成されるキャリア数で表されています。光捕捉構造によって吸収が促進されている箇所では明るい領域が現れ、活性層内に生成が集中している様子がはっきりと確認できます。この生成分布データこそが電気モデルへと渡される実際の入力値となるのです。

電流密度電圧特性。最後に、デバイス性能を迅速に評価する上で最も重要な点として、jvを開いて??に示されている電流密度電圧特性を確認してください。この例では、セルの開放電圧（\(V_{oc}\)）は1 V強、短絡電流密度（\(J_{sc}\)）は約−150 A m−2となっています。負の符号は、光生成によって生じる電流が順方向バイアス電圧の方向とは逆に流れるという符号規約を反映したものです。なお、OghmaNanoでは全ての物理量がSI単位で表記される点にご留意ください。

6. 結果の詳細な検証

今回実行したシミュレーションは実質的に二つの部分から構成されています。すなわち、FDTDによる光学計算とドリフト拡散法による電気的特性計算です。OghmaNanoはこれら両段階における詳細なスナップショットを保存しており、それらを順次確認することでデバイス内部で実際にどのような現象が起きているかを把握できます。まずは光学計算の結果を確認し、次に電気的特性計算の結果を見ていきます。

FDTD（光学）スナップショットについて説明します。[Output]タブからfdtd_snapshotsをダブルクリックすると、??に示すスナップショットビューアーが開きます。これはシミュレーション実行中に得られたFDTD電場データの完全なダンプ内容となっています。このビューアーを使えば、シミュレーションステップごとに結果を確認できます。青色の＋ボタンをクリックしてpower_density.csvを選択した後、スクロールバーを動かすと、電場がデバイス内部へどのように伝播し構造物によって散乱されていくかを観察できます。もし生の電場成分そのものを確認したい場合は、パワー密度の代わりにEx、Ey、またはEzを追加すればよいのです。

電気的特性のスナップショット機能について。電気的特性の計算結果は別途保存される。Outputタブからsnapshotsディレクトリをダブルクリックする（fdtd_snapshotsではない）。さらに＋ボタンを使用してEc.csv、Ev.csv、Fn.csv、Fp.csvを追加することで、伝導帯端、価電子帯端、および電子・正孔の擬フェルミ準位をプロットできる。これらは??に示されている。なお今回はスライダーが時間ではなく印加電圧を変化させるため（ウィンドウタイトルには現在のバイアス値、例えば0.1 Vが表示される）、短絡状態から開路状態へと電圧を変化させた際のバンド構造の応答を観察することが可能となる。

キャリア密度。これと密接に関連する表示としては、Q_nfree.csvおよびQ_pfree.csvを追加することで得られるものがあり、これらはそれぞれ自由電子密度および自由正孔密度を表す（単位は\(\mathrm{m^{-3}}\)）。その結果が??に示されている。同様に、スライダーを電圧値に沿って動かすことで、印加バイアスが変化するにつれてデバイス内で電子および正孔の分布がどのように再配分されるかを確認することができる。

7. 光捕捉構造の効果

セル上部にある三角形状の層は、光閉じ込め構造となっており（同時に反射防止層の役割も果たす）。この構造がどのような効果をもたらすかを正確に測定するため、シミュレーションにおいてこの層を無効化し、元のデバイスの結果と比較することが可能である。

基準となる電流電圧特性を保存する。何も変更する前に、結果ファイル jv.csv のコピーを作成し、オペレーティングシステムのファイルマネージャーを使用してそのコピーに jv_backup.csv という名前を付ける。これにより、元のデバイスの電流電圧特性が保持され、後で比較する際に利用できる。

構造体を無効化します。デバイスの上部にある三角形構造体を右クリックし、??に示すように「オブジェクト編集」を選択します。これにより??に示すオブジェクトエディタが開きます。ツールバーにある有効/無効ボタンをクリックし、「オブジェクト無効」と表示される状態にします。そうすると三角形構造体がデバイスから消え、??に示す平坦な積層構造のみが残ります。

実際に我々が無効化したのは何か？このようなテクスチャ層は太陽電池にとって二つの利点をもたらす。第一に、急峻な平坦面を傾斜した表面に置き換えることで表面反射が抑制され、より多くの光が素子内部へと入り込めるようになる。第二に、透過した光が散乱・方向転換され、吸収層内をより長い経路を辿ることで吸収の可能性が高まる。これらの効果により全体的な光吸収量、ひいては短絡電流が増加するのだ。実際のデバイスにおいては、例えばシリコン太陽電池の表面に形成されるミクロンスケールのピラミッド構造がその典型例であり、この効果は広く認知されているが実験的に分離・定量化するのは困難な場合も多い。シミュレーションを用いればこの構造を有効・無効に切り替え、その影響を直接的に把握することが可能となる。

再度シミュレーションを実行し比較してみましょう。三角形構造を無効化した状態でシミュレーションを再実行すると、光捕捉層のないデバイス用の新たなjv.csvが得られます。この新しいjv.csvを、保存済みのjv_backup.csv（光捕捉層を備えた完全なデバイス）と同一の座標軸上にプロットすると??のような結果が得られるはずです。その差はわずかですが明確に認識できます。つまり、光捕捉構造がある場合の短絡電流の絶対値は約−162 A m−2であり、構造がない場合の約−155 A m−2よりも大きくなっています。これにより、このテクスチャ構造が吸収される光量を増加させ、結果としてセルが供給できる電流も増大することが実証されます。

8. 電気的問題の次元削減

多次元シミュレーションを実施する必要は必ずしもない。多くの多層構造を持つ太陽電池においては、主な変化はデバイスの厚み方向に生じるため、十分な解像度を持つ一次元モデルでも重要な物理現象を低い計算コストで再現可能である。多次元モデルが最も有用となるのは、デバイスの幾何形状やキャリア輸送が横方向にも変化する場合である。

この例では、テクスチャー加工された光閉じ込め表面の影響により、光学的な解析は実質的に二次元問題となる。光は側方へ散乱されるため、電磁場はx方向およびデバイス厚み方向の両方で変化する。一方、電気的解析は異なる様相を呈する。すなわち、電荷輸送は依然として各層を通じた移動が支配的であるため、一次元のドリフト拡散モデルでも十分な精度が得られる可能性がある。

この違いは極めて重要である。なぜなら、2次元FDTD法と2次元ドリフト拡散法を組み合わせたシミュレーションは、単純化されたモデルに比べて計算コストがはるかに高くなるからだ。そこで我々は、表面テクスチャの影響が顕著な光学計算のみを2次元で行い、電気計算については目的とする特性を損なうことなく1次元モデルへ簡略化できるかどうかを検証することにする。

??に示す電気解析用リボンを開き、[Electrical mesh]をクリックします。そうすると??に示すメッシュエディターが表示されます。このデバイスは、層構造に沿ってy方向に40点、横方向のx方向に10点のメッシュで離散化されています。より細かいy方向のメッシュにより、バンド曲がり、キャリア勾配、再結合現象、およびデバイス内の電流の流れが正確に解析可能となります。

電気モデルを一次元化するには、メッシュエディター上部にある「X Dimension」をクリックします。これにより横方向のxメッシュが削除され、厚み方向のyメッシュのみが残ります。その様子は??に示されています。光学的シミュレーションは依然として二次元のままですが、ドリフト拡散計算はより高速な一次元電気モデルとして実行されるようになります。

シミュレーションを再実行してください。このモデルでは、2次元光学シミュレーションが1次元電気シミュレーションと連動しています。FDTDソルバーは依然として光捕捉構造によって生じる横方向の散乱を解析しますが、ドリフト拡散ソルバーはキャリア輸送をデバイスの厚み方向のみで取り扱います。

これまで通り電流密度電圧特性を確認してください。次に snapshots ディレクトリを開き、自由電子および正孔の密度をプロットします。その結果は ?? と類似するはずです。完全な電気モデルによって得られる2次元面分布の代わりに、今回は層構造に沿った1次元プロファイルとしてキャリア密度が表示されます。

これはしばしば効率的なモデリング手法となります。すなわち、幾何形状上の必要性から2Dモデルを採用しつつ、物理現象が主に厚み方向にのみ作用する場合には1Dモデルを用いるのです。この近似法は、横方向のキャリア輸送が重要でない場合に適しており、横方向に均一な薄膜太陽電池においてはしばしば該当します。

9. 光学的問題の次元削減

電気的な問題の次元を低減させたので、今度は光学的な問題についても同様の処理を行い、シミュレーション全体をどれだけコンパクトにできるかを確認しましょう。??に示す光学リボンに移動し、光学メッシュをクリックします。これにより??に示す光学メッシュエディターが開きます。電気メッシュの場合と同様に、X次元ボタンをクリックしてX方向の離散化を解除します。これにより光学的問題は??に示すように一次元で解かれることになり、光学計算の規模が大幅に縮小されるのです。

光学メッシュエディターを開いている間に、右側にある波長メッシュにも目を向けてみましょう。現在のシミュレーションでは、波長範囲300～800ナノメートルの中間付近、およそ550ナノメートルという単一の波長点のみで計算が行われていることがわかります（これは先ほどのFDTD計算で用いた緑色光と同じ波長です）。つまり、この例では太陽光スペクトル全体をシミュレートしているわけではなく、単一色の光に対してのみデバイスの特性を評価しているのです。もし太陽光スペクトル全体に対する応答を調べたい場合は、より多くの波長点を追加する必要があります（スペクトルを適切に解析するためには少なくとも20点程度は必要でしょう）。それは十分に可能かつ妥当な方法ですが、計算時間が大幅に増加することになるため、ここではデモの実行速度を維持するために単一波長のみを採用しています。

さて、再びシミュレーションを実行してください。これで1次元光学シミュレーションが1次元電気シミュレーションと結合され、最も軽量なモデルが構築されます。??に示すように、スナップショットビューアーから1次元FDTDの結果を確認できます。パルスがデバイス内を伝播する際、光パワー密度はチュートリアル開始時に見られた2次元の場分布ではなく、単一の1次元線状プロファイルとなります。

10. より高速な代替手法：伝達行列法

チュートリアルを完了するため、これまで説明した計算負荷の非常に高いFDTD法から、はるかに計算コストの低い伝達行列法（TMM）へと切り替えます。多くの太陽電池シミュレーションツールにおいて伝達行列法は純粋に一次元モデルと見なされており、平坦な平面構造にしか適用できないと考えられています。しかしOghmaNanoの構造上、任意の方向・任意の次元において伝達行列法を適用可能であり、多くの場合においてFDTD法の柔軟な代替手段となり得るのです。

光学メニューに移動し、転送行列ボタンをクリックしてください。これにより??で示される光学シミュレーションエディタが開きます。続いて「光学シミュレーション実行」をクリックすると、そこに示されているような光子分布が得られます。OghmaNanoではこの際、2次元デバイス全体に対して転送行列法を適用し、構造を複数の独立した垂直スライスに分割して解析しています。重要な点として、FDTDとは異なり隣接するスライス間での光の相互結合は生じません。転送行列法はマクスウェル方程式のより単純な近似手法であるため、光があるスライスから隣接するスライスへ横方向に散乱することはあり得ないのです。この近似の代償として、電場の時間ステップ計算が不要なため、FDTDよりもはるかに高速にシミュレーションが実行されるのです。

転送行列法は計算コストが非常に低いため、単一波長だけでなく多数の波長に対してもシミュレーションを実行することが可能となる。光学メッシュエディターに戻り、波長メッシュを約20点に設定した上で光学シミュレーションを再実行してください。実行中には、??に示すようにターミナル上で順次各波長が処理されていく様子が確認できる。各光学スライスはそれぞれ固有の波長（ここでは約312.5 nmから787.5 nmまで）で解析され、現在処理中の波長を示すためにテキストが色分けされている。

計算が完了したら、転送行列ウィンドウ??に戻り、「吸収された光子分布」タブを選択してください。ウィンドウ下部にあるスライダーを操作することで波長を順次切り替え、各色の光がデバイス内のどこで吸収されているかを確認できます。例ではスライダーが587.5 nmに設定されており、ペロブスカイトデバイス内の対応する吸収分布が表示されています。この方法により、太陽光スペクトルの各成分がセル内のどの深さで吸収されるかを迅速かつ直感的に把握することが可能です。

11. 光生成モデルの選択

デフォルト設定では、本シミュレーションではFDTDソルバーが光生成モデルとして使用されており、これが電気計算に用いられる生成分布を供給する。より高速な転送行列法を電気シミュレーションのための光生成モデルとして利用したい場合は、生成モデルを直接変更することが可能だ。??に示すリボンを開き、JVエディターをクリックする。JVエディター内の「電荷キャリア生成」項目において、FDTDから転送行列法へと電荷キャリア生成モデルを変更する（??参照）。これ以降は、光生成分布がFDTDではなく転送行列法によって算出されるようになる。

まとめ。本チュートリアルを通じて、我々は単一のデバイスに対してさまざまなモデリング手法を適用してきた。まず、光捕捉構造を有するペロブスカイト太陽電池に対し、2次元FDTD法による光学シミュレーションと2次元ドリフト拡散法による電気的特性解析を連携させた。次に、光学および電気的問題の次元を2次元から1次元へと低減させ、計算時間や結果の変化を観察した後、再び2次元モデルへと戻した。最後に、FDTD法に代わりはるかに高速な転送行列法を用いることで、全く別のアプローチで光学問題の簡略化を実現した。