Secs1d package

From Octave
Revision as of 12:25, 2 September 2018 by Siko1056 (talk | contribs) (Rename "Octave-Forge" to "Octave Forge" (https://lists.gnu.org/archive/html/octave-maintainers/2018-08/msg00138.html).)
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to navigation Jump to search
  • Example, a PIN diode

% physical constants and parameters
 secs1d_physical_constants;
 secs1d_silicon_material_properties;
 
 % geometry
 L  = 10e-6;          % [m] 
 xm = L/2;
 
 Nelements = 1000;
 x         = linspace (0, L, Nelements+1)';
 sinodes   = [1:length(x)];
 
 % dielectric constant (silicon)
 er = esir * ones (Nelements, 1);
 
 % doping profile [m^{-3}]
 Na = 1e23 * (x <= xm);
 Nd = 1e23 * (x > xm);
 
 % avoid zero doping
 D  = Nd - Na;  
  
 % initial guess for n, p, V, phin, phip
 V_p = -1;
 V_n =  0;
 
 Fp = V_p * (x <= xm);
 Fn = Fp;
 
 p = abs (D) / 2 .* (1 + sqrt (1 + 4 * (ni./abs(D)) .^2)) .* (x <= xm) + ...
     ni^2 ./ (abs (D) / 2 .* (1 + sqrt (1 + 4 * (ni ./ abs (D)) .^2))) .* (x > xm);
 
 n = abs (D) / 2 .* (1 + sqrt (1 + 4 * (ni ./ abs (D)) .^ 2)) .* (x > xm) + ...
     ni ^ 2 ./ (abs (D) / 2 .* (1 + sqrt (1 + 4 * (ni ./ abs (D)) .^2))) .* (x <= xm);
 
 V = Fn + Vth * log (n / ni);
 
 % scaling factors
 xbar = L;                       % [m]
 nbar = norm(D, 'inf');          % [m^{-3}]
 Vbar = Vth;                     % [V]
 mubar = max (u0n, u0p);         % [m^2 V^{-1} s^{-1}]
 tbar = xbar^2 / (mubar * Vbar); % [s]
 Rbar = nbar / tbar;             % [m^{-3} s^{-1}]
 Ebar = Vbar / xbar;             % [V m^{-1}]
 Jbar = q * mubar * nbar * Ebar; % [A m^{-2}]
 CAubar = Rbar / nbar^3;         % [m^6 s^{-1}]
 abar = 1/xbar;                  % [m^{-1}]
 
 % scaling procedure
 l2 = e0 * Vbar / (q * nbar * xbar^2);
 theta = ni / nbar;
 
 xin = x / xbar;
 Din = D / nbar;
 Nain = Na / nbar;
 Ndin = Nd / nbar;
 pin = p / nbar;
 nin = n / nbar;
 Vin = V / Vbar;
 Fnin = Vin - log (nin);
 Fpin = Vin + log (pin);
 
 tnin = tn / tbar;
 tpin = tp / tbar;
 
 u0nin = u0n / mubar;
 uminnin = uminn / mubar;
 vsatnin = vsatn / (mubar * Ebar);
 
 u0pin = u0p / mubar;
 uminpin = uminp / mubar;
 vsatpin = vsatp / (mubar * Ebar);
 
 Nrefnin = Nrefn / nbar;
 Nrefpin = Nrefp / nbar;
 
 Cnin     = Cn / CAubar;
 Cpin     = Cp / CAubar;
 
 anin     = an / abar;
 apin     = ap / abar;
 Ecritnin = Ecritn / Ebar;
 Ecritpin = Ecritp / Ebar;
 
 % tolerances for convergence checks
 toll  = 1e-3;
 maxit = 1000;
 ptoll = 1e-12;
 pmaxit = 1000;
 
 % solve the problem using the full DD model
 [nout, pout, Vout, Fnout, Fpout, Jnout, Jpout, it, res] = ...
       secs1d_dd_gummel_map (xin, Din, Nain, Ndin, pin, nin, Vin, Fnin, Fpin, ...
                             l2, er, u0nin, uminnin, vsatnin, betan, Nrefnin, ...
 	                       u0pin, uminpin, vsatpin, betap, Nrefpin, theta, ...
 		               tnin, tpin, Cnin, Cpin, anin, apin, ...
 		               Ecritnin, Ecritpin, toll, maxit, ptoll, pmaxit); 
 
 % Descaling procedure
 n    = nout*nbar;
 p    = pout*nbar;
 V    = Vout*Vbar;
 Fn   = V - Vth*log(n/ni);
 Fp   = V + Vth*log(p/ni);
 dV   = diff(V);
 dx   = diff(x);
 E    = -dV./dx;
 
 % band structure
 Efn  = -Fn;
 Efp  = -Fp;
 Ec   = Vth*log(Nc./n)+Efn;
 Ev   = -Vth*log(Nv./p)+Efp;
 
 plot (x, Efn, x, Efp, x, Ec, x, Ev)
 legend ('Efn', 'Efp', 'Ec', 'Ev')
 axis tight