% Berechnung der Elektronendichte n0 bei n-dot. in Abhängigkeit v.d.
% Temperatur

clear all;
close all;

ND1 = 10^14; % aus Aufgabenstellung
ND2 = 10^16; %  "        "
Wg  = 1.08;  % [eV]
WV  = 0;
z   = 125;  % 1/Zoomfaktor
T   = logspace(1,3,500);

% Berechnen der Zustandsdichten
NC  = 2.81E19 * (T/300).^(3/2);
NV  = 1.83E19 * (T/300).^(3/2);
k   = 8.62E-5;
kT  = k*T;

% Festlegen des Nullpunktes
WC  = WV+Wg;
WC_ = WC./kT;
Wg_ = Wg./kT;

% Donatorniveau berechnen
WD  = WC-0.024-kT*log(2);
WD_ = WD./kT;

% Eigenleitungsdichte S.93 Gl. 2.34
ni  = sqrt(NC.*NV).*exp(-Wg_/2);

% Seite 113 Gl. 2.60
nob1= 2*ND1./( 1+ sqrt(1+4* ND1./NC .* exp(WC_-WD_)));
nob2= 2*ND2./( 1+ sqrt(1+4* ND2./NC .* exp(WC_-WD_)));
no1 = ni + nob1;
no2 = ni + nob2;
no  = [no1; no2];
x   = 1000*1./T;

% Ergebnisse plotten
figure()
semilogy(x,no1,'--',x,no2,'-.','LineWidth',3)
hold on
grid on
legend('N_D = 10^{14}', 'N_D = 10^{16}', 'location','south')
ylabel('Ladungstraeger Anzahl')
xlabel('1/T')
figure()
semilogy(x(length(x)-z:length(x)),no1(length(x)-z:length(x)),'--',x(length(x)-z:length(x)),no2(length(x)-z:length(x)),'-.','LineWidth',3)
ylabel('Ladungstraeger Anzahl')
xlabel('1/T')
legend('N_D = 10^{14}', 'N_D = 10^{16}', 'location','northeast')
grid on