Modelowanie rozproszonych stanów pułapkowych w nieuporządkowanych półprzewodnikach.

Niniejszy poradnik stanowi ogólne wprowadzenie do modelowania stanów pułapkowych w nieuporządkowanych półprzewodnikach. PM6:Y6 wykorzystano jako przykład, ponieważ jest to dobrze znany organiczny materiał fotowoltaiczny, ale te same zasady modelowania można stosować w szerszym zakresie do organicznych półprzewodników, amorficznego krzemu, półprzewodników tlenkowych, perowskitów i innych materiałów cienkowarstwowych, w których zlokalizowane stany tworzą rozproszoną gęstość stanów wewnątrz przerwy energetycznej.

1. Wprowadzenie: dlaczego rozkłady stanów pułapek są istotne w nieuporządkowanych półprzewodnikach.

Stanów pułapkowych odgrywają kluczową rolę w zachowaniu wielu nieuporządkowanych półprzewodników, w tym półprzewodników organicznych, amorficznego krzemu, półprzewodników tlenkowych, perowskitów i innych materiałów fotowoltaicznych w postaci cienkich warstw. Różnicę między uporządkowanymi i nieuporządkowanymi półprzewodnikami można zobaczyć na ??. W wysoce uporządkowanych półprzewodnikach, takich jak krzem krystaliczny o wyjątkowo wysokiej czystości, często opisywany jako materiał o czystości rzędu „sześć dziewiątek” lub „dziewięć dziewiątek”, występuje stosunkowo niewiele stanów pułapkowych w przerwie energetycznej. Te izolowane stany pułapkowe są zilustrowane na obszarze a na ??. W takich układach rekombinację zachodzącą przez te defekty można często opisać za pomocą standardowego równania rekombinacji Shockleya–Reada–Halla (SRH):

\[ R_{SRH}=\frac{np-n_i^2} {\tau_p(n+n_1)+\tau_n(p+p_1)} \]

Szczegółowe, krok po kroku przedstawione wyprowadzenie matematyczne znajduje się w Wyprowadzenie równania Shockleya–Reada–Halla dotyczącego rekombinacji.

To równanie zakłada, że rekombinacja zachodzi poprzez pojedynczy, dyskretny stan pułapki w środkowej części przerwy energetycznej, co często stanowi odpowiednie przybliżenie dla stosunkowo uporządkowanych półprzewodników/urządzeń, takich jak diody p–n z krzemu i diody p–n z arsenku galu. Jednak wraz ze wzrostem stopnia nieuporządkowania, coraz mniej przydatne staje się postrzeganie pułapek jako izolowanych defektów, a bardziej przydatne jest opisywanie ich jako całego rozkładu stanów wewnątrz przerwy energetycznej. Ten rozkład stanów jest powszechnie określany jako gęstość stanów pułapki lub gęstość stanów (DOS). Prawa strona ?? ilustruje tę sytuację dla nieuporządkowanego półprzewodnika, takiego jak półprzewodnik organiczny lub amorficzny krzem. Obszar b reprezentuje rozkład stanów pułapki zawierających elektrony, znajdujących się blisko pasma przewodnictwa, obszar d reprezentuje rozkład stanów pułapki zawierających dziury, znajdujących się blisko pasma walencyjnego, a obszar c reprezentuje rozkład głębszych stanów pułapki w środkowej części przerwy energetycznej.

W uporządkowanych półprzewodnikach standardowe równanie Shockleya–Reada–Halla (SRH) jest często wystarczające, ponieważ stany pułapkowe są zazwyczaj rzadkie i można je często przybliżyć jako izolowane poziomy defektów. W nieuporządkowanych półprzewodnikach to przybliżenie przestaje być trafne: materiał zawiera szeroki rozkład zlokalizowanych stanów rozciągających się wzdłuż przerwy energetycznej, więc założenie, że pułapki stanowią pojedynczy punktowy defekt, jest niewystarczające. Równanie SRH dla pojedynczego poziomu traktuje pułapkę jako jeden dyskretny poziom energii, dlatego nie może opisać, w jaki sposób nośniki wypełniają szeroki rozkład zlokalizowanych stanów. W nieuporządkowanym półprzewodniku to wypełnienie określa równowagę między ruchomym ładunkiem a ładunkiem uwięzionym w pułapkach, a zatem ustala zależność między całkowitą gęstością nośników a położeniem quasi-poziomu Fermiego. Jeśli ta zależność jest nieprawidłowa, model podaje również nieprawidłową szybkość rekombinacji, efektywną ruchliwość i charakterystykę J–V. Aby poprawnie modelować nieuporządkowane urządzenia, należy wyraźnie rozwiązać równania dotyczące wypełnienia rozkładu pułapek, wykorzystując równania szybkości pułapkowania, uwalniania z pułapki i rekombinacji, które są sprzężone w sposób samospójny z równaniami dryfu-dyfuzji i równaniem Poissona. Otrzymana populacja uwięzionych nośników nie jest małą poprawką: może dominować w procesach magazynowania ładunku, transportu, rekombinacji, ruchliwości, prądu w stanie ciemnym, napięcia w obwodzie otwartym oraz kształcie fotowoltaicznej krzywej J–V.

2. Zjawisko wychwytywania i uwalniania nośników ładunku w nieuporządkowanych półprzewodnikach, opisywane przez model Shockleya–Reada–Halla.

Wielu czytelników zapoznaje się z procesem rekombinacji Shockleya–Reada–Halla (SRH) po raz pierwszy dzięki przedstawionemu powyżej zwartemu równaniu opisującemu stan ustalony. Jednak to zwarte równanie jest jedynie uproszczonym wynikiem wyprowadzonym z bardziej rozbudowanego formalizmu równań szybkości SRH; poszczególne kroki matematyczne przedstawiono w Wyprowadzenie rekombinacji Shockleya–Reada–Halla. W istocie teoria SRH opisuje wypełnianie i opróżnianie stanu pułapkowego poprzez procesy wychwytywania i uwalniania nośników, co zilustrowano w ??. Elektrony mogą być wychwytywane do pułapki, a następnie emitowane z niej pod wpływem temperatury i wracać do pasma przewodnictwa, podczas gdy dziury również mogą być wychwytywane i emitowane. Gdy zarówno elektron, jak i dziura znajdują się w tym samym stanie pułapkowym, następuje rekombinacja.

Zatem podstawowa formalna struktura SRH jest zasadniczo modelem opartym na równaniach kinetycznych, opisującym wychwytywanie, uwalnianie i rekombinację nośników ładunku:

\[ \frac{dn_t}{dt}=r_{ec}-r_{ee}-r_{hc}+r_{he} \]

W tym kontekście \(r_{ec}\) oznacza szybkość wychwytu elektronów przez pułapkę, \(r_{ee}\) – szybkość ucieczki elektronów z powrotem do pasma przewodnictwa, \(r_{hc}\) – szybkość wychwytu dziur, a \(r_{he}\) – szybkość ucieczki dziur. Jedną z kluczowych zasad teorii SRH jest to, że szybkości ucieczki (\(r_{ee}\) i \(r_{he}\)) silnie zależą od głębokości energetycznej pułapki, natomiast szybkości wychwytu (\(r_{ec}\) i \(r_{hc}\)) zależą głównie od tego, czy pułapka jest zajęta. Płytkie pułapki uwalniają nośniki stosunkowo łatwo, natomiast głębokie pułapki mogą magazynować ładunek przez długi czas.

OghmaNano rozwiązuje pełne formalne równanie dotyczące pułapkowania nośników, zamiast stosować upraszczające założenia, które zazwyczaj wykorzystuje się do wyprowadzenia zwartego równania SRH opisującego rekombinację w stanie ustalonym. W praktyce oznacza to, że dla każdego poziomu pułapkowania w obrębie gęstości stanów rozwiązywane jest oddzielne równanie pułapkowania. Dzięki temu symulator może modelować zgromadzone ładunki, wypełnianie pułapek, rekombinację wspomaganą przez pułapki, elektrostatyczne zakrzywienie pasm energetycznych oraz rozkład populacji pułapek w silnie nieuporządkowanych półprzewodnikach, takich jak półprzewodniki organiczne i amorficzny krzem.

W tym samouczku omówimy, w jaki sposób zaimplementowano i rozwiązano problem rozproszonych stanów pułapkowych w programie OghmaNano, oraz w jaki sposób procesy zatrzymywania, uwalniania i akumulacji ładunku wpływają na działanie nieuporządkowanych urządzeń półprzewodnikowych.

3. Przykład: tworzenie modelu symulacji półprzewodnika o nieuporządkowanej strukturze PM6:Y6.

Uruchom program OghmaNano i kliknij „Nowa symulacja”. Otworzy się biblioteka dostępnych typów urządzeń, jak pokazano na ilustracji ??. Kliknij dwukrotnie ikonę „Organiczne ogniwa słoneczne”, aby otworzyć bibliotekę przykładów OPV. Z listy dostępnych przykładów PM6:Y6 wybierz PM6:Y6_E10_0hrs, jak pokazano na ilustracji ??. W tym samouczku można również użyć ogólnej symulacji dryfu i dyfuzji; jednak ten przykład PM6:Y6 został już skonfigurowany z rozkładem stanów pułapkowych i strukturą pasmową, które odpowiadają wynikom przedstawionym tutaj. Po wyświetleniu odpowiedniego komunikatu zapisz symulację w folderze, do którego masz uprawnienia zapisu.

Po utworzeniu symulacji otwiera się główne okno programu OghmaNano, jak pokazano na ilustracji ??. Urządzenie PM6:Y6 jest przedstawione jako warstwowa struktura ogniwa słonecznego. Warstwa aktywna to nieuporządkowany obszar półprzewodnikowy, w którym rozproszone stany pułapkowe zostaną zmodyfikowane.

4. Ustawienie symulacji tak, aby działała w ciemności.

Zanim rozpoczniemy symulację, musimy najpierw skonfigurować urządzenie, aby działało w ciemności. Otwórz pasek optyczny, który znajduje się w ??, i ustaw natężenie światła na 0,0 jednostek słonecznych. Zapewni to, że generowane krzywe JV będą odpowiadały pomiarom JV w ciemności, a nie działaniu ogniwa fotowoltaicznego w warunkach oświetlenia. Symulacje w ciemności są często szczególnie przydatne podczas badania stanów pułapkowych, ponieważ ujawniają w sposób bardziej wyraźny wstrzykiwanie nośników, wypełnianie pułapek, rekombinację i efekty elektrostatyczne niż symulacje w warunkach oświetlenia, gdzie fotoprąd może maskować wiele z tych procesów.

5. Określenie gęstości stanów w pułapce.

W głównym oknie otwórz edytor parametrów elektrycznych dla aktywnej warstwy PM6:Y6, korzystając z zakładki dotyczącej struktury urządzenia. Spowoduje to otwarcie edytora przedstawionego na rysunku ??. Gdy włączona jest opcja „Dynamiczne pułapki SRH”, OghmaNano oblicza rozkład stanów pułapek w przestrzeni energii, który przypomina rozkład zilustrowany na rysunku ??.

Domyślnie zakłada się, że rozkład pułapek odpowiada wykładniczej gęstości stanów. Dla obszaru o niskiej energii LUMO gęstość stanów wyraża się wzorem:

\[ \rho_\mathrm{LUMO}(E)=N_{t,n}\exp\left(\frac{E-E_C}{E_{U,n}}\right), \]

natomiast w przypadku ogona po stronie HOMO zapisuje się go jako:

\[ \rho_\mathrm{HOMO}(E)=N_{t,p}\exp\left(\frac{E_V-E}{E_{U,p}}\right). \]

W powyższym tekście \(N_{t,n}\) i \(N_{t,p}\) oznaczają odpowiednio współczynniki preeksponencjalne dla gęstości pułapek elektronów i dziur, natomiast \(E_{U,n}\) i \(E_{U,p}\) oznaczają odpowiadające im energie w ogonie Urbacha. W interfejsie użytkownika przedstawionym na ilustracji ?? odpowiadają one odpowiednio parametrom: gęstość pułapek elektronów, gęstość pułapek dziur, nachylenie ogona elektronów i nachylenie ogona dziur w sekcji „Pułapki SRH w stanie nieustalonym” edytora parametrów elektrycznych.

Edytor zawiera również przekroje przechwytujące, które kontrolują szybkość procesów wychwytywania nośników ładunku, o których mowa wcześniej. Parametry „Swobodny elektron do elektronu uwięzionego” i „Swobodna dziura do elektronu uwięzionego” regulują procesy przechwytywania elektronów i dziur (\(r_{ec}\) i \(r_{hc}\)), natomiast parametry „Swobodny elektron do dziury uwięzionej” i „Swobodna dziura do dziury uwięzionej” regulują odpowiednie procesy uwalniania (\(r_{ee}\) i \(r_{he}\)). Razem parametry te określają, jak szybko nośniki ładunku są przechwytywane i uwalniane z rozkładu pułapek.

Na koniec, parametr „Liczba pułapek” kontroluje liczbę dyskretyzowanych poziomów pułapek, które są wykorzystywane do przybliżenia analitycznej gęstości stanów. Odpowiadają one poziomom pułapek, które są przedstawione jako poziome linie w kolorze czerwonym i niebieskim na rysunku ??. Sama analityczna gęstość stanów jest nadal definiowana przez funkcje \(\rho_\mathrm{LUMO}(E)\) i \(\rho_\mathrm{HOMO}(E)\), natomiast liczba pułapek określa, jak dokładnie rozkład jest dyskretyzowany numerycznie.

Jeśli symulacja zostanie uruchomiona w tej konfiguracji, OghmaNano przyjmie wykładniczą zależność gęstości stanów, odpowiadającą powyższym równaniom. OghmaNano umożliwia jednak również znacznie większą kontrolę nad dokładnym rozkładem stanów pułapkowych. Jeśli zostanie kliknięty przełącznik „Wykładniczy”, tryb przełącza się na „Złożony”, jak pokazano na ilustracji ??. Pozwala to na zdefiniowanie dowolnych kształtów matematycznych gęstości stanów. Następnie można otworzyć edytor złożonej gęstości stanów, klikając przycisk „Edytuj”, który pojawia się po prawej stronie edytora.

6. Zdefiniuj rozkład płytkiej eksponencjalnej pułapki.

Kliknij „Edytuj”, aby otworzyć zaawansowany edytor gęstości stanów. Edytor jest wyświetlany na ??. Górna tabela definiuje rozkład pułapek związanych z poziomem LUMO, a dolna tabela definiuje rozkład pułapek związanych z poziomem HOMO. Każdy wiersz zawiera funkcję matematyczną, przełącznik aktywacji oraz trzy parametry numeryczne oznaczone jako a, b i c. Wykres po lewej stronie edytora pokazuje wynikową gęstość stanów, generowaną na podstawie funkcji matematycznych wprowadzonych do tabel. Niebieska krzywa analityczna reprezentuje rozkład pułapek związanych z poziomem LUMO lub pasmem przewodnictwa (\(E_C\)), natomiast czerwona krzywa analityczna reprezentuje rozkład pułapek związanych z poziomem HOMO lub pasmem walencyjnym (\(E_V\)). Te gładkie krzywe odpowiadają dokładnej analitycznej gęstości stanów, zdefiniowanej przez równania wprowadzone do edytora.

Niebiesko-zielone, zakrzywione schody symbolizują dyskretyzowane poziomy pułapek, które są faktycznie wykorzystywane wewnętrznie przez symulator. Te poziomy schodów bezpośrednio odpowiadają parametrowi „Liczba pułapek”, zdefiniowanemu wcześniej w edytorze parametrów elektrycznych. Zwiększenie liczby poziomów pułapek prowadzi do dokładniejszego przybliżenia analitycznej gęstości stanów, a tym samym poprawia dokładność numeryczną symulacji. W praktyce symulacje wykorzystujące około pięciu poziomów pułapek zazwyczaj dają bardzo dobre wyniki dla wielu symulacji urządzeń w stanie ustalonym. Jednak symulacje obejmujące długie okresy przejściowe, takie jak obliczenia czasu przelotu, mogą wymagać większej liczby poziomów pułapek, czasami nawet do około dwudziestu. Ponieważ stany pułapek są rozwiązywane w każdym punkcie siatki w obrębie urządzenia, zwiększenie liczby poziomów pułapek bardzo szybko zwiększa ogólną wielkość matematyczną symulacji, dlatego w praktyce rzadko wymagane są wartości znacznie większe niż dwadzieścia.

Co do zasady, dowolne równanie można wprowadzić do edytora równań, korzystając z parametrów zawartych w powyższej tabeli, a następnie można wykorzystać to do określenia gęstości stanów.

7. Wpływ szerokości obszaru „pułapki” na charakterystykę prądowo-napięciową w ciemności oraz na gęstość ładunku.

Przeprowadzimy teraz dwie symulacje i porównamy wpływ zmiany energetycznej szerokości rozkładu pułapek. W zaawansowanym edytorze gęstości stanów upewnij się, że parametr b jest ustawiony na 0,025 eV zarówno dla rozkładu LUMO, jak i HOMO, zgodnie z ilustracją w ??. Określa to przypadek płytkiego opadania, w którym gęstość stanów pułapek maleje stosunkowo szybko w miarę oddalania się od krawędzi pasma, a zatem w przerwie energetycznej występuje tylko ograniczona liczba głębokich stanów pułapek.

Uruchamianie symulacji i tworzenie kopii zapasowej do wykorzystania w przyszłości: Wróć do głównego okna i uruchom symulację, klikając niebieski przycisk. Uruchom symulację. przycisk lub poprzez naciśnięcie F9Po zakończeniu obliczeń otwórz folder symulacji za pomocą Eksploratora Windows. Plik jv.csv zawiera charakterystyczną, ciemną krzywą JV, natomiast charge.csv zawiera informacje o całkowitym ładunku zgromadzonym w urządzeniu, przedstawione w zależności od napięcia. Utwórz kopie zapasowe tych plików za pomocą Eksploratora Windows i zmień nazwy skopiowanych plików na: jv.bak i charge.bakTe pliki kopii zapasowych będą zawierały wyniki symulacji „płaskiego ogona”, dzięki czemu w późniejszym terminie będzie można je porównać z szerszym rozkładem pułapek.

Następnie ponownie otwórz zaawansowany edytor DOS i zwiększ wartość parametru b do 0,100 eV zarówno dla rozkładu LUMO, jak i HOMO, zgodnie z ilustracją ??. Powoduje to poszerzenie rozkładu stanów pułapkowych, dzięki czemu stany pułapkowe rozciągają się dalej w obszar przerwy energetycznej. Uruchom ponownie symulację. Po drugim uruchomieniu w katalogu symulacji powinny znajdować się zarówno oryginalne pliki kopii zapasowych (jv.bak i charge.bak), jak i nowo wygenerowane pliki wynikowe (jv.csv i charge.csv).

Po przeprowadzeniu obu symulacji można użyć programu Microsoft Excel lub poniższych skryptów MATLAB i Python, aby wygenerować wykresy przedstawione w ?? i ??. W tych porównaniach pliki .bak odpowiadają symulacji z płytkim obszarem opadania, natomiast nowe pliki .csv odpowiadają symulacji z szerszym obszarem opadania.

Zwiększenie energii „ogona Urbacha” powoduje poszerzenie zakresu gęstości stanów pułapkowych, dzięki czemu znacznie większa liczba pułapek rozciąga się głęboko w obszar przerwy energetycznej. To bezpośrednio wpływa zarówno na prąd rekombinacji, jak i na ilość ładunku zgromadzonego w pułapkach wewnątrz urządzenia. Na wykresie porównawczym charakterystyk prądowo-napięciowych w ciemności, przedstawionym w ??, szerszy rozkład pułapek zwiększa prąd w średnim zakresie napięć, ponieważ nośniki ładunku mogą uzyskać dostęp do znacznie większej liczby stanów pułapkowych o różnej energii, które biorą udział w rekombinacji wspomaganej przez pułapki i w transporcie ładunku.

Efekt ten jest również wyraźnie widoczny w porównaniu rozkładu gęstości ładunku, przedstawionym na ??. Szerszy rozkład pułapek powoduje, że w całym urządzeniu gromadzona jest większa ilość ładunku, ponieważ teraz dostępnych jest znacznie więcej stanów energetycznych, w których mogą znajdować się nośniki ładunku. Nagromadzony ładunek działa jako ładunek przestrzenny wewnątrz urządzenia, co zmienia potencjał elektrostatyczny i zakrzywienie pasm energetycznych.

Z punktu widzenia fizyki, to jest zasadnicza różnica między modelowaniem pojedynczego centrum rekombinacji SRH a modelowaniem pełnego, rozproszonego rozkładu stanów pułapkowych. W nieuporządkowanym półprzewodniku, pułapki nie są izolowanymi defektami. Tworzą one całe spektrum energetyczne, które jednocześnie kontroluje wychwytywanie i uwalnianie nośników, rekombinację, magazynowanie ładunku i transport. Dlatego też stany pułapkowe muszą być uwzględnione w modelach urządzeń półprzewodnikowych o nieuporządkowanej strukturze.

clear;
close all;
clc;

% ============================================
% Settings
% ============================================

log_jv = true;

clear;
close all;
clc;

% ============================================
% Settings
% ============================================

log_jv = true;    % true = semilogy
                  % false = linear plot

% ============================================
% Load data
% ============================================

jv         = load('jv.csv');
jv_bak     = load('jv.bak');

charge     = load('charge.csv');
charge_bak = load('charge.bak');

% ============================================
% Figure 1 : JV
% ============================================

figure;

if log_jv == true

    semilogy(jv(:,1), abs(jv(:,2)), ...
        'LineWidth', 2);

    hold on;

    semilogy(jv_bak(:,1), abs(jv_bak(:,2)), '--', ...
        'LineWidth', 2);

else

    plot(jv(:,1), jv(:,2), ...
        'LineWidth', 2);

    hold on;

    plot(jv_bak(:,1), jv_bak(:,2), '--', ...
        'LineWidth', 2);

end

grid on;
box on;

xlabel('Voltage (V)');
ylabel('Current density (A m^{-2})');

title('JV comparison');

legend( ...
    'JV', ...
    'JV backup', ...
    'Location', 'northwest' ...
);

% ============================================
% Figure 2 : Charge density
% ============================================

figure;

semilogy(charge(:,1), abs(charge(:,2)), ...
    'LineWidth', 2);

hold on;

semilogy(charge_bak(:,1), abs(charge_bak(:,2)), '--', ...
    'LineWidth', 2);

grid on;
box on;

xlabel('Voltage (V)');
ylabel('Charge density (m^{-3})');

title('Charge density comparison');

legend( ...
    'Charge', ...
    'Charge backup', ...
    'Location', 'northwest' ...
);

import numpy as np
import matplotlib.pyplot as plt

# ============================================
# Settings
# ============================================

log_jv = True

import numpy as np
import matplotlib.pyplot as plt

# ============================================
# Settings
# ============================================

log_jv = True    # True = semilogy
                 # False = linear plot

# ============================================
# Load data
# ============================================

jv = np.loadtxt("jv.csv")
jv_bak = np.loadtxt("jv.bak")

charge = np.loadtxt("charge.csv")
charge_bak = np.loadtxt("charge.bak")

# ============================================
# Figure 1 : JV
# ============================================

plt.figure()

if log_jv:
    plt.semilogy(jv[:, 0], np.abs(jv[:, 1]), linewidth=2)
    plt.semilogy(jv_bak[:, 0], np.abs(jv_bak[:, 1]), "--", linewidth=2)
else:
    plt.plot(jv[:, 0], jv[:, 1], linewidth=2)
    plt.plot(jv_bak[:, 0], jv_bak[:, 1], "--", linewidth=2)

plt.grid(True)

plt.xlabel("Voltage (V)")
plt.ylabel("Current density (A m$^{-2}$)")
plt.title("JV comparison")
plt.legend(["JV", "JV backup"], loc="upper left")

# ============================================
# Figure 2 : Charge density
# ============================================

plt.figure()

plt.semilogy(charge[:, 0], np.abs(charge[:, 1]), linewidth=2)
plt.semilogy(charge_bak[:, 0], np.abs(charge_bak[:, 1]), "--", linewidth=2)

plt.grid(True)

plt.xlabel("Voltage (V)")
plt.ylabel("Charge density (m$^{-3}$)")
plt.title("Charge density comparison")
plt.legend(["Charge", "Charge backup"], loc="upper left")

plt.show()