Kinetic and Macroscopic Models for Semiconductors

(1)

Kinetic and Macroscopic Models for Semiconductors

Ansgar J¨ungel

Vienna University of Technology, Austria

www.jungel.at.vu

Ansgar J¨ungel (TU Wien) Kinetic Semiconductor Models www.jungel.at.vu 1 / 165

(2)

Literature

Main reference

A. J¨ungel. Transport Equations for Semiconductors. Springer, 2009.

Physics of semiconductors:

K. Brennan. The Physics of Semiconductors. Cambridge, 1999.

M. Lundstrom. Fundamentals of Carrier Transport. Cambridge, 2000.

Kinetic semiconductor models:

P. Markowich, C. Ringhofer, and C. Schmeiser. Semiconductor Equations. Vienna, 1990.

Macroscopic semiconductor models:

P. Degond. Mathematical modelling of microelectronics semiconductor devices. Providence, 2000.

F. Brezzi, L. Marini, S. Micheletti, P. Pietra, R. Sacco, and S. Wang.

Discretization of semiconductor device problems. Amsterdam, 2005.

(4)

From motherboard to transistor

CopyrightMoxfyre aten.wikipedia

Processor

Copyright Ioan Sameli

Electric circuit

Copyright FDominec

Transistor

Copyright Daniel Ryde

(5)

History of Intel processors

1971

4004

108 KHz, 2250 transistors,

feature size: 10µm (1µm= 10⁻⁶m)

1982

80286

12 MHz, 134,000 transistors, feature size: 1.5µm

1996

Pentium 2

66 MHz, 7,500,000 transistors, feature size: 0.35µm

2011

Core i7

3.5 GHz, 2,600,000,000 transistors, feature size: 0.032µm= 32nm

(6)

Channel lengths 2000–2016

(7)

Challenges in semiconductor device design

Future processors (2012):

Number of transistors>2,000,000,000 Transistor feature size 22 nm

Highly-integrated circuits:

power density >100 W/cm² Some key problems:

Decreasing power supply → noise effects

Increasing frequencies → multi-scale problems

Increasing design variety → need of fast and accurate simulations Increasing power density → parasitic effects (heating, hot spots)

(8)

What are semiconductors?

Non-conducting at temperature T = 0 K, conducting at T >0 (heat, light etc.)

Modern definition: energy gap of order of a few eV Basis material: Silicon, Germanium, GaAs etc.

Doping of the basis material with other atoms, gives higher conductivity

Modeled by doping concentrationC(x)

(9)

How does a semiconductor transistor work?

Elektronen

Source Gate Drain

Bulk

70 Nanometer

MOSFET = Metal-Oxide Semiconductor Field-Effect Transistor Source and drain contact: electrons flow from source to drain Gate contact: applied voltage controls electron flow

Advantage: small gate voltage controls large electron current Used as an amplifier or switch

(10)

Objectives

Introduce into the basics of semiconductor physics

Modeling of electron transport through semiconductors by kinetic equations

Modeling of macroscopic electron transport (numerically cheaper than kinetic models)

Modeling of quantum transport and quantum diffusion effects Numerical approximation of macroscopic models (finite-element and finite-difference methods)

(11)

Overview

5 Summary

(12)

What do we need from physics?

Newton’s law:

x(t): spatial variable, v(t): velocity, m: mass,F: force

˙

x=v, mv˙ =F, t >0 Schr¨odinger equation:

Schr¨odinger equation for wave functionψ:

i~∂_tψ=−~²

2m∆ψ−qV(x,t)ψ in R³ where~=h/2π,q: elementary charge,V(x,t): potential Interpretation of wave function:

particle density: n=|ψ|², current density: J=−~

mIm(ψ∇ψ) Stationary equation: ψ(x,t) =e^−iEt/~φ(x) gives eigenvalue problem

−~²

2m∆φ−qV(x)φ=Eφ, E : energy

(13)

Classical particle transport

Given ensemble of M particles of massmmoving in a vacuum

Trajectory (x(t),v(t))∈R³×R³ computed from Newton’s equations

˙

x =v, mv˙ =F, t >0, x(0) =x0, v(0) =v0

Force: F =∇V(x,t),V(x,t): electric potential

M 1: use statistical description with probability density f(x,v,t) Theorem (Liouville)

Let x˙ =X(x,v),v˙ =V(x,v). If

∂X

∂x +∂V

∂v = 0 then f(x(t),v(t),t) =f_I(x₀,v₀), t >0

→ Assumption satisfied if F =F(x,t)

(14)

Vlasov equation

Differentiation off(x(t),v(t),t) =f_I(x₀,v₀) givesVlasov equation:

0 = d

dtf(x(t),v(t),t) =∂tf +x˙ · ∇xf +v˙ · ∇vf

=∂tf +v· ∇xf + 1

m∇xV(x,t)· ∇vf Moments of f(x,v,t):

Particle density: n(x,t) = Z

R³

f(x,v,t)dv Current density: J(x,t) =

Z

R³

vf(x,v,t)dv Energy density: (ne)(x,t) =

Z

R³

m

2|v|²f(x,v,t)dv Electrons are quantum mechanical objects: quantum description needed!

(15)

Electrons in a semiconductor

Semiconductor = ions (nuclei + core electrons) and valence electrons State of ion-electron system described by wave functionψ

Schr¨odinger eigenvalue problem:

−~²

2m∆ψ−qV_L(x)ψ=Eψ, x ∈R³ VL=Vei +Veff: periodic lattice potential

V_ei: electron-ion Coulomb interactions

V_eff: effective electron-electron interactions (Hartree-Fock approx.) Goal: exploit periodicity of lattice potential

(16)

Electrons in a semiconductor

Schr¨odinger eigenvalue problem:

−~²

2m∆ψ−qVL(x)ψ=Eψ, x∈R³ Theorem (Bloch)

Schr¨odinger eigenvalue problem in R³ can be reduced to Schr¨odinger problem on lattice cell, indexed by k ∈B (B: dual cell or Brillouin zone)

−~²

2m∆ψ−qV_L(x)ψ=Eψ, x ∈cell

For each k, there exists sequence(E, ψ) = (E_n(k), ψ_n,k), n ∈N ψ_n,k(x) =e^ik·xu_n,k(x), where u_n,k periodic on lattice

En(k) is real, periodic, symmetric on Brillouin zone En(k) =n-th energy band

energy gap = allE^∗ for which there is no k with En(k) =E^∗

(17)

Energy bands

Silicon Gallium Arsenide

(18)

Parabolic band approximation

Taylor expansion aroundk = 0 if E(0) = 0:

E(k)≈E(0) +∇kE(0)·k+1

2k^>d²E dk²(0)k

= 1

2k^>d²E dk²(0)k Diagonalization:

1

~² d²E

dk²(0) =





1/m^∗₁ 0 0

0 1/m^∗₂ 0

0 0 1/m^∗₃





isotropic

=





1/m^∗ 0 0

0 1/m^∗ 0

0 0 1/m^∗





Parabolic band approximation

E(k) = ~² 2m^∗|k|²

(19)

Semi-classical picture

i~∂tψ=−~²

2m∆ψ−q(VL(x) +V(x))ψ where V_L: lattice potential, V: external potential

Theorem (Semi-classical equations of motion)

~x˙ =~vn(k) =∇kEn(k), ~k˙ =q∇xV Momentum operator: Pψ_n,k = (~/i)∇ψ_n,k

Mean velocity: vn=hPi/m= (~/im)R

ψ_n,k∇ψn,kdx Motivation of the formulas:

Insert ψ_n,k(x) =e^ik^·xu_n,k(x) in Schr¨odinger equation⇒ first eq.

ψ_n,k(x) =e^ik·x ⇒ Pψ_n,k =~kψ_n,k: ~k = crystal momentum =p Newton’s law: ~k˙ = ˙p =F =q∇xV gives second equation

(20)

Effective mass

Semi-classical equations of motion:

~x˙ =~v_n(k) =∇kE_n(k), ~k˙ =q∇xV Definition of effective massm^∗:

p =m^∗v_n Consequence:

˙

p=m^∗ ∂

∂tvn= m^∗

~

∂

∂t∇kEn= m^∗

~ d²En

dk² k˙ = m^∗

~² d²En

dk² p˙ Effective mass equation:

m^∗ =~² d²En

dk² −1

(21)

Overview

5 Summary

(22)

Semi-classical kinetic equations

Semi-classical equations:

~x˙ =∇kE(k), ~k˙ =q∇xV(x), p =m^∗v Liouville’s theorem: If(and this is true)

∂

∂x∇kE(k) + ∂

∂kq∇xV(x) = 0 then f(x(t),k(t),t) =f_I(x₀,k₀) Semi-classical Vlasov equation:

0 = d

dtf(x,k,t) =∂_tf+ ˙x·∇xf+ ˙k·∇kf =∂_tf+v(k)·∇xf+q

~∇xV·∇kf Include collisions: assume thatdf/dt =Q(f)

Semi-classical Boltzmann equation

∂_tf +v(k)· ∇xf +q

~∇xV · ∇kf =Q(f)

(23)

Poisson equation

Electric force given by E =E_ext+E_mean

Electric force between electrons given by Coulomb field:

Ec(x,y) =− q 4πεs

x−y

|x−y|³ Mean-field approximation of electric field:

Emean(x,t) = Z

R³

n(y,t)Ec(x,y)dy ⇒ divEmean=−q εs

n External electric field generated by doping atoms:

E_ext(x,t) = q 4πεs

Z

R³

C(y) x−y

|x−y|³dy ⇒ divE_ext= q εs

C(x) Since curlE = 0, there exists potential V such thatE =−∇V Poisson equation

εs∆V =−εsdiv(Emean+Eext) =q(n−C(x))

(24)

Holes

Hole = vacant orbital in valence band

Interpret hole as defect electron with positive charge

Current flow = electron flow in conduction band and hole flow in valence band

Electron density n(x,t), hole density p(x,t)

(25)

Holes

- E(k) 6

k d

t

?

valence band conduction band

energy emission -

- E(k) 6

k t

6 d

valence band conduction band

energy absorption

Recombination: conduction electron recombines with valence hole Generation: creation of conduction electron and valence hole Shockley-Read-Hall model:

R(n,p) = n²_i −np

τp(n+n_d) +τn(p+p_d), ni : intrinsic density

(26)

Boltzmann distribution function

∂tf +v(k)· ∇xf + q

~∇xV · ∇kf =Q(f), v(k) =∇kE(k)/~

Definition of distribution function:

f(x,k,t) = number of occupied states in dx dk in conduction band total number of states in dx dk in conduction band Quantum state has phase-space volume (2π)³ (integrate

k ∈B ∼(−π, π)³)

Quantum state density (take into accountelectron spin):

N^∗(x,k)dx dk= 2

(2π)³dx dk = 1 4π³dx dk Total number of electrons in volumedx dk:

dn=f(x,k,t)N^∗(x,k)dx dk =f(x,k,t)dx dk 4π³ Electron density:

n(x,t) = Z

B

dn= Z

B

f(x,k,t) dk 4π³

(27)

Collision models

Probability that electron changes statek⁰ to k is proportional to occupation prob. f(x,k⁰,t)× non-occupation prob.(1−f(x,k,t)) Collisions between two electrons in states k andk⁰:

(Q(f))(x,k,t) = (Probabilityk⁰ →k)−(Probability k →k⁰)

= Z

B

s(x,k⁰,k)f⁰(1−f)−s(x,k,k⁰)f(1−f⁰) dk⁰ wheref⁰ =f(x,k⁰,t),s(x,k⁰,k): scattering rate

Important collision processes:

Electron-phonon scattering Ionized impurity scattering Electron-electron scattering

(28)

Scattering rates

Electron-phonon scattering:

Collisions of electrons with vibrations of crystal lattice (phonons) Phonon emission: E(k⁰)−E(k) =~ω = phonon energy

Phonon absorption: E(k⁰)−E(k) =−~ω

Phonon occupation number: N = 1/(exp(~ω/k_BT)−1) General scattering rate:

s(x,k,k⁰) =σ (1 +N)δ(E⁰−E+~ω) +Nδ(E⁰−E −~ω) whereδ: delta distribution,E⁰ =E(k⁰)

If phonon scattering elastic,~ω ≈0: s(x,k,k⁰) =σ(x,k,k⁰)δ(E⁰−E) (Q_el(f))(x,k,t) =

Z

B

σ(x,k,k⁰)δ(E⁰−E)(f⁰−f)dk⁰ Elastic collisions conserve mass and energy:

Z

B

Q_el(f)dk = Z

B

E(k)Q_el(f)dk = 0

(29)

Scattering rates

Ionized impurity scattering:

Collisions of electrons with ionized doping atoms: elastic scattering Collision operator

(Q(f))(x,k,t) = Z

B

σ(x,k,k⁰)δ(E⁰−E)(f⁰−f)dk⁰ Electron-electron scattering:

Electrons in states k⁰ andk₁⁰ collide and scatter to statesk andk₁ Elastic collisions: s(x,k,k⁰,k₁,k₁⁰) =σ(x)δ(E⁰+E₁⁰−E−E₁) Collision operator:

(Q(f))(x,k,t) = Z

B³

s(x,k,k⁰,k₁,k₁⁰)

× f⁰f₁⁰(1−f)(1−f₁)−ff₁(1−f⁰)(1−f₁⁰)

dk⁰dk₁dk₁⁰ Mass and energy conservation: R

BQ(f)dk =R

BE(k)Q(f)fdk = 0

(30)

Summary

Electron motion in semi-classical approximation:

∂tf +v(k)· ∇xf +q

~∇xV · ∇kf =Q(f), x ∈R³, k ∈B B: Brillouin zone coming from crystal structure

k: pseudo-wave vector, p =~k: crystal momentum Mean velocity: v(k) =∇kE(k)/~

Energy band E(k); parabolic band approximation:

E(k) =~²|k|²/2m^∗

Electric potentialV computed from Poisson equation ε_s∆V =q(n−C(x)), C(x) : doping profile Electron density:

n(x,t) = Z

B

f(x,k,t) dk 4π³

(31)

Overview

5 Summary

(32)

Why macroscopic models?

∂_tf +v(k)· ∇xf +q

~∇xV · ∇kf =Q(f), x∈R³, k ∈B ⊂R³ Semi-classical Boltzmann equation is (3+3)-dimensional: numerical simulations extremely time-consuming

Often only macroscopic physical variables (particle density n, velocity u, energy ne) are of interest

(n,nu,ne) = Z

R³

f(k)(1,k,¹₂|k|²) dk 4π³

Derive evolution equations by integrating Boltzmann equation over k ∈B

Depending on semiconductor application, derive various models→ leads to model hierarchy

(33)

Scaling of Boltzmann equation

∂_tf +v(k)· ∇xf +q

~∇xV · ∇kf =Q(f), ε_s∆V =q(n−C(x)) Introduce reference values for

length λ time τ

mean free path λc =uτ =λ velocity u=p

k_BT_L/m^∗, wave vector k₀ =m^∗u/~ potential U_T =k_BT_L/q Scaled Boltzmann equation:

∂tf +v(k)· ∇xf +∇xV · ∇kf =Q(f) Scaled Poisson equation:

λ²_D∆V =n−C(x), λ²_D = εsU_T qλ²k₀

Objective: derive macroscopic equations by averaging over k ∈B

(34)

Moment method

Boltzmann equation with parabolic band:

∂_tf +k· ∇xf +∇xV · ∇kf =Q(f), λ²_D∆V =n−C(x) Integrate over k ∈B:

∂t

Z

B

f dk 4π³

| {z }

=n(x,t)

+divx

Z

B

kf dk 4π³

| {z }

=−Jn(x,t)

+∇^xV · Z

B∇kf dk 4π³

| {z }

=0

= Z

B

Q(f) dk 4π³

| {z }

=0

→ Mass balance equation: ∂tn−divJ_n= 0 Multiply by k and integrate by parts:

∂t

Z

B

kf dk 4π³

| {z }

=−J_n(x,t)

+divx

Z

B

k⊗kf dk 4π³

| {z }

=P

−∇xV · Z

B

f dk 4π³

| {z }

=n

= Z

B

kQ(f) dk 4π³

| {z }

=−W

→ Momentum balance equation: ∂_tJ_n−divP+∇V ·J_n =W

(35)

Moment method: closure problem

Mass balance equation:

∂_tn−divJ_n= 0 Momentum balance equation:

∂_tJ_n−divP +∇V ·J_n=W, P = Z

B

k⊗kf dk 4π³

Energy balance equation (assuming energy conservation):

∂t

Z

B

|k|² 2 f dk

4π³

| {z }

=(ne)(x,t)

+divZ

B

k|k|² 2

dk 4π³

| {z }

=R

−∇V· Z

B

kf dk 4π³

| {z }

=−Jn

= Z

B

|k|²

2 Q(f) dk 4π³

| {z }

=0

→ ∂t(ne) +divR+∇V ·J_n= 0

Closure problem: P and R cannot be expressed in terms of n,Jn,ne

(36)

Solution of closure problem

Scaling of Boltzmann equation:

Collision timeτ_c =τ /α: hydrodynamic scaling (long time scale) α∂_tf+α v·∇xf+∇xV·∇kf

=Q(f), α= λ_c

λ = Knudsen number Collision timeτc =τ /α²: diffusion scaling (very long time scale) α²∂tf+α v·∇^xf+∇^xV·∇kf

=Q(f), Equilibrium distribution (Maxwellian):

Kinetic entropy: S(f) =−R

Bf(logf −1 +E(k))dk Given f, solve constrained maximization problem:

max n

S(g) : Z

B

κ(k)g dk 4π³ =

Z

B

κ(k)f dk

4π³, κ(k) = 1,k,1 2|k|²o

→ formal solution: Maxwellian M[f] = exp(κ(k)·λ(x))

(37)

Model hierarchy

Model hierarchy depends on . . . diffusive or hydrodynamic scaling number of moments or weight functions Hydrodynamic models:

Weight functions 1, k: isothermal hydrodynamic equations for electron density n and current densityJ_n

Weight functions 1, k, ¹₂|k|²: full hydrodynamic equations forn,Jn, and energy density ne

Diffusive models:

Weight function 1: drift-diffusion equations forn

Weight functions 1, ¹₂|k|²: energy-transport equations for n andne

(38)

Model hierarchy

Diffusive models Hydrodynamic models Drift-diffusion

equations

Isothermal hydrodynamic equations Energy-transport

equations

Full hydrodynamic equations Fourth-order moment

equations

Extended hydrodynamic equations Higher-order moment

equations Stratton 1962 Van Roosbroeck 1950

Anile 1995 Grasser et al. 2001

Blotekjaer 1970

Higher-order hydrodynamic equations Struchtrup 1999 A.J./Krause/Pietra

2007

1

2

3

# Variables

4

5

13

(39)

Overview

5 Summary

(40)

Derivation

α²∂tf_α+α v(k)· ∇xf_α+∇xV · ∇kf_α

=Q(f_α) Simplifications: parabolic bandE(k) = ¹₂|k|² (k ∈R³), relaxation-time operatorQ(f) = (nM−f)/τ

Maxwellian: M(k) = (2π)^−3/2exp(−¹₂|k|²), R

R³M(k)dk = 1 Electron density: n_α(x,t) =R

R³f_α(x,k,t)dk/4π³ Moment equation: integrate Boltzmann equation overk

α∂t

Z

R³

fα

dk

4π³ + divx

Z

R³

kfα

dk 4π³ = 1

ατ Z

R³

(nM−fα)dk 4π³ Derivation in three steps

Step 1: limitα→0 in Boltzmann equation⇒ Q(f) = 0

⇒ f = limα→0fα =nM

(41)

Derivation

Step 2:

Chapman-Enskog expansion fα=nαM+αgα in Boltzmann equation:

α∂tfα+ k· ∇x(nαM) +∇xV · ∇k(nαM) +α k· ∇xg_α+∇xV · ∇kg_α

=α⁻¹Q(n_αM) +Q(g_α) =Q(g_α) Limitα→0 (g = limα→0g_α): Q(g) = (nM−g)/τ

Q(g) =k· ∇x(nM) +∇xV · ∇k(nM) =k·(∇xn−n∇xV)M

⇒ g =−τk·(∇xn−n∇xV)M +nM, M(k) = (2π)^−3/2e^−|k^|²^/2 Step 3:

Insert Chapman-Enskog expansion in Boltzmann equation:

∂_t Z

R³

f_α dk 4π³+1

αdiv_x Z

R³

kn_αM dk 4π³

| {z }

=0

+div_x Z

R³

kg_α dk 4π³ = 1

α²τ Z

R³

Q(f_α) dk 4π³

| {z }

=0

(42)

Derivation

∂t

Z

R³

(nαM+αgα) dk

4π³ + divx

Z

R³

kgα

dk 4π³ = 0 Limitα→0:

∂_t Z

R³

nM dk

4π³ + div_x Z

R³

kg dk 4π³ = 0 Define current density J_n=−R

R³kgdk/4π³, insert expression for g =−τk·(∇xn−n∇xV)M+nM:

∂_tn−divJ_n= 0, J_n=τ Z

R³

k⊗kM dk 4π³

| {z }

=Id

(∇xn−n∇xV)

Theorem (Drift-diffusion equations) The formal limitα→0 gives

∂tn−divJn= 0, Jn=τ(∇n−n∇V)

(43)

Bipolar drift-diffusion equations

∂tn−divJ_n=−R(n,p), J_n=τ(∇n−n∇V)

∂_tp+divJ_p =−R(n,p), J_n=−τ(∇p+p∇V) λ²_D∆V =n−p−C(x)

Hole density modeled by drift-diffusion equations Shockley-Read-Hall recombination-generation term:

R(n,p) = np−n²_i

τp(n+n_d) +τn(p+p_d) with physical parameter n_i,τ_n,τ_p,n_d,p_d

Auger recombination-generation term (high carrier densities):

R(n,p) = (Cnn+Cpp)(np−n²_i) with physical parameter Cn and Cp

Equilibrium state: np=n²_i = intrinsic density

(44)

Drift-diffusion equations: summary

∂tn−divJn= 0, Jn=∇n−n∇V, λ²_D∆V =n−C(x) Variables: electron densityn, electric potential V

n∇V: drift current, ∇n: diffusion current First proposed by van Roosbroeck 1950

Rigorous derivation from Boltzmann equation: Poupaud 1992 (linear), Ben Abdallah/Tayeb 2004 (1D Poisson coupling), Masmoudi/Tayeb 2007 (multi-dimensional)

Existence analysis: Mock 1972, Gajewski/Gr¨oger 1986

Numerical solution: Scharfetter/Gummel 1964, Brezzi et al. 1987 + well established, used in industrial semiconductor codes

+ well understood analytically and numerically + stable mixed finite-element schemes available

− satisfactory results only for lengths>1µm

− no carrier heating (thermal effects)

(45)

Drift-diffusion models: extensions

Assumption: Generalization:

¹ Equilibrium given by Maxwellian M Use Fermi-Dirac distribution ² Parabolic-band approximation Use non-parabolic bands ³ Scaling valid for low fields only Devise high-field models ⁴ One moment (particle density) used Use more moments Comparison drift-diffusion and Monte-Carlo:

From: Grasser et al. 2005

(46)

Summary

Standard

drift-diffusion model High-field drift-diffusion model Drift-diffusion using

Fermi-Dirac statistics

High-density drift-diffusion model Drift-diffusion models

diffusion approximations Non-parabolic band

drift-diffusion model parabolic band

(47)

Overview

5 Summary

(48)

Semi-classical Boltzmann equation

∂tf +v(k)· ∇xf +q

~∇xV · ∇kf =Q(f) Collision operator: Q(f) =Qel(f) +Qee(f) +Qin(f)

Q_el(f) = Z

B

σ_el(k,k⁰)δ(E⁰−E)(f⁰−f)dk⁰ Qee(f) =

Z

B³

σee(,k,k⁰,k1,k₁⁰)δ(E⁰+E₁⁰−E−E1)

× f⁰f₁⁰(1−f)(1−f₁)−ff₁(1−f⁰)(1−f₁⁰)

dk⁰dk₁dk₁⁰ Qin(f) = inelastic collisions (unspecified)

Scaling: α=p λ_el/λ_in

α²∂tf +α v(k)· ∇xf +∇xV · ∇kf

=Q_el(f) +αQee(f) +α²Qin(f)

(49)

Properties of elastic collision operator

Q_el(f) = Z

B

σ_el(k,k⁰)δ(E⁰−E)(f⁰−f)dk⁰, σ(k,k⁰) symmetric Proposition

Conservation properties: R

BQ_el(f)dk =R

BQ_el(f)E(k)dk = 0 for all f Symmetry: −Qel is symmetric and nonnegative

Kernel N(Qel)= all functions F(x,E(k),t) Proof:

Conservation and symmetry: use symmetry of σ(k,k⁰) andδ(E⁰−E) Nonnegativity: show that

Z

B

Qel(f)fdk = 1 2

Z

B²

σel(k,k⁰)δ(E⁰−E)(f⁰−f)²dk⁰dk ≥0 Kernel: Q_el(f) = 0⇒ δ(E⁰−E)(f⁰−f)²= 0 ⇒ f(k⁰) =f(k) if E(k⁰) =E(k) ⇒ f constant on energy surface {k :E(k) =ε}

(50)

Properties of elastic collision operator

Q_el(f) = Z

B

σ_el(k,k⁰)δ(E⁰−E)(f⁰−f)dk⁰, σ(k,k⁰) symmetric Proposition

Equation Q_el(f) =h solvable iffR

Bh(k)δ(E(k)−ε)dk = 0 for all ε Proof:

Fredholm alternative: Q_el symmetric⇒ Q_el(f) =h solvable iff h∈N(Q_el)^⊥

Let Q_el(f) =h be solvable and leth∈N(Q_el)^⊥,f =F(E)∈N(Q_el):

0 = Z

B

hfdk = Z

B

h(k) Z

R

F(ε)δ(E(k)−ε)dεdk

= Z

R

Z

B

h(k)δ(E(k)−ε)dk F(ε)dε ⇒ Z

B

h(k)δ(E(k)−ε)dk = 0 Conversely, show similarly that h∈N(Q_el)^⊥

(51)

Properties of electron-electron collision operator

Q_ee(f) = Z

B³

σee(k,k⁰,k₁,k₁⁰)δ(E⁰+E₁⁰ −E −E₁)

× f⁰f₁⁰(1−f)(1−f1)−ff1(1−f⁰)(1−f₁⁰)

dk⁰dk1dk₁⁰ Proposition

Let σee be symmetric:

BQee(f)dk =R

BQee(f)E(k)dk = 0 ∀f Kernel N(Qee) = Fermi-Dirac distributions F(k),

F(k) = 1/(1 + exp((E(k)−µ)/T)) for arbitraryµ, T Proof: Show that

Z

B

Q_ee(f)gdk =− Z

B⁴

σ_eeδ(E +E₁−E⁰−E₁⁰)(g⁰+g₁⁰ −g −g₁)

× f⁰f₁⁰(1−f)(1−f1)−ff1(1−f⁰)(1−f₁⁰) dk⁴ Conservation: takeg = 1 andg =E. Kernel: more difficult

(52)

Properties of electron-electron collision operator

Qee(f) = Z

B³

σee(k,k⁰,k1,k₁⁰)δ(E⁰+E₁⁰ −E −E1)

× f⁰f₁⁰(1−f)(1−f1)−ff1(1−f⁰)(1−f₁⁰)

dk⁰dk1dk₁⁰ Averaged collision operator:

S(ε) = Z

B

Q_ee(F)δ(E−ε)dk

Proposition

Let σee be symmetric:

RS(ε)dε=R

RS(ε)εdε= 0

If S(ε) = 0for all εthen F = 1/(1 + exp((E −µ)/T))Fermi-Dirac Proof: similar as above

(53)

Derivation: general strategy

α²∂_tf_α+α v(k)· ∇xf_α+∇xV · ∇kf_α

=Q_el(f_α) +αQ_ee(f_α) +α²Q_in(f_α) Set hgi=R

Bg(k)dk/4π³

Moment equations for momentshf_αi andhEf_αi:

α²∂_thE^jf_αi+αdiv_xhE^jvf_αi −α∇xV · h∇kE^jf_αi

=hE^jQ_el(f_α)i+αhE^jQ_ee(f_α)i+α²hE^jQ_in(f_α)i

=α²hE^jQ_in(f_α)i, j = 0,1 Strategy of derivation:

Step 1: formal limit α→0 in Boltzmann equation Step 2: Chapman-Enskog expansion f_α=F +αg_α Step 3: formal limit α→0 in moment equations

References: Ben Abdallah/Degond 1996, Degond/Levermore/Schmeiser 2004

(54)

Step 1

α²∂tfα+α v(k)· ∇^xfα+∇^xV · ∇kfα

=Q_el(fα) +αQee(fα) +α²Qin(fα) Step 1: α→0 in Boltzmann equation⇒ Q_el(f) = 0, wheref = limα→0fα

⇒ f(x,k,t) =F(x,E(k),t) Step 2:

Chapman-Enskog expansion f_α=F +αgα in Boltzmann equation:

α∂_tf_α+ v(k)· ∇xF +∇xV · ∇kF +α v(k)· ∇xg_α+∇xV · ∇kg_α

=Q_el(g_α) +Q_ee(f_α) +αQ_in(f_α) Formal limit α→0 gives

Qel(g) =v(k)· ∇xF +∇xV · ∇kF −Qee(F) Operator equation solvable iff

Z

B

v(k)· ∇^xF +∇^xV · ∇kF −Qee(F)

δ(E−ε)dk = 0 ∀ε

(55)

Step 2

Solvability condition for operator equation:

Z

B

v(k)· ∇xF +∇xV · ∇kF −Qee(F)

δ(E−ε)dk = 0 ∀ε Since ∇kF =∂_EF∇kE,v =∇kE andH⁰ =δ (H: Heaviside fct.)

Z

B

v(k)· ∇xF +∇xV · ∇kF

δ(E−ε)dk

= ∇^xF +∂_EF∇^xV (ε)·

Z

B∇kEδ(E−ε)dk

= ∇^xF +∂_EF∇^xV (ε)·

Z

B∇kH(E−ε)dk = 0 Solvability condition becomes

Z

B

Qee(F)

δ(E−ε)dk = 0

⇒ F_µ,T = Fermi Dirac

(56)

Step 3

Operator equation becomes (withF_µ,T = 1/(1 + exp((E(k)−µ)/T))) Q_el(g) =v(k)· ∇^xF_µ,T +∇^xV · ∇kF_µ,T −Qee(F_µ,T)

| {z }

=0

=F_µ,T(1−F_µ,T)v(k)·

∇x

µ

T −∇xV

T −E∇x

1 T

Step 3: limitα→0 in the moment equations Set hgi=R

Bg(k)dk/4π³. Moment equations for j = 0,1:

∂thE^jfαi+α⁻¹hE^j(v· ∇xFµ,T +∇xV · ∇kFµ,T)i

| {z }

=hE^jQel(g)i=0

+hE^j(v· ∇xg_α+∇xV · ∇kg_α)i=hE^jQ_in(f_α)i Limitα→0:

∂_thE^jFi+hE^j(v· ∇xg+∇xV · ∇kg)i

| {z }

=divxhE^jvgi−∇xV·h∇kE^jgi

=hE^jQ_in(F)i

(57)

Step 3: balance equations

Moment equation forj = 0:

(assume mass conservation for inelastic scattering)

∂t hFi

|{z}

=n

+divxhvgi

|{z}

=−J₀

−∇xV · h∇kgi

| {z }

=0

=hQin(F)i

| {z }

=0

Moment equation forj = 1:

∂_thEFi

| {z }

=ne

+div_xhEvgi

| {z }

=−J1

−∇xV · h∇kEgi

| {z }

=−J0

=hE^jQ_in(F)i

| {z }

=W

Particle current density J₀ =−hvgi Energy current densityJ₁ =−hEvgi Energy relaxation term W =hEQin(F)i Balance equations

∂tn−divJ0 = 0, ∂t(ne)−divJ1+∇V ·J0=W

(58)

Step 3: current densities

J₀=−hvgi, J₁=−hEvgi where g is solution of

Qel(g) =Fµ,T(1−Fµ,T)v(k)·

∇µ T −∇V

T −E∇1 T

Let d0 be solution ofQel(d0) =−Fµ,T(1−Fµ,T)v(k). Then g =−d0·

∇µ T −∇V

T −E∇1 T

+F1, F1 ∈N(Q_el) Insert into expressions for current densities:

J₀ =D₀₀

∇µ T −∇V

T

−D₀₁∇1 T J1 =D10

∇µ T −∇V

T

−D11∇1 T Diffusion coefficients:

Dij =hE^i+jv⊗d0i= Z

B

E^i+jv⊗d0

dk 4π³

(59)

Summary

Energy-transport equations

∂tn−divJ0 = 0, ∂t(ne)−divJ1+∇V ·J0=W, x ∈R³, t >0 J₀ =D₀₀

∇µ T −∇V

T

−D₀₁∇1

T, J₁=D₁₀

∇µ T −∇V

T

−D₁₁∇1 T Electron and energy densities:

n(µ,T) = Z

B

F_µ,T dk

4π³, ne(µ,T) = Z

B

E(k)F_µ,T dk 4π³ Diffusion coefficients:

Dij = Z

B

E^i+jv⊗d0

dk

4π³, d0 solves Qel(d0) =−Fµ,T(1−Fµ,T)v Energy-relaxation term:

W(µ,T) = Z

B

E(k)Qin(F_µ,T) dk 4π³

(60)

References

∂tn−divJ0 = 0, ∂t(ne)−divJ1+∇V ·J0=W, x ∈R³, t >0 J₀ =D₀₀

∇µ T −∇V

T

−D₀₁∇1

T, J₁=D₁₀

∇µ T −∇V

T

−D₁₁∇1 T First energy-transport model: Stratton 1962 (Rudan/Gnudi/Quade 1993)

Derivation from Boltzmann equation: Ben Abdallah/Degond 1996 Existence results:

Heuristic temperature model: Allegretto/Xie 1994

Uniformly positive definite diffusion matrix: Degond/G´enieys/A.J. 1997 Close-to-equilibrium solutions: Chen/Hsiao/Li 2005

Numerical approximations:

Mixed finite volumes: Bosisio/Sacco/Saleri/Gatti 1998

Mixed finite elements: Marrocco/Montarnal 1996, Degond/A.J./Pietra 2000, Holst/A.J./Pietra 2003-2004, Gadau/A.J. 2008

(61)

Relation to nonequilibrium thermodynamics

∂_tn−divJ₀ = 0, ∂_t(ne)−divJ₁+∇V ·J₀=W, x ∈R³, t >0 J0 =D00

∇µ T −∇V

T

−D01∇1

T, J1=D10

∇µ T −∇V

T

−D11∇1 T Balance equations = conservation laws of mass and energy (if no forces)

Thermodynamic forces:

X₀ =∇(µ/T)− ∇V/T, X₁ =−∇(1/T) Thermodynamic fluxes:

J0 =D00X0+D01X1, J1 =D10X0+D11X1

Density variables n,ne Entropy variables µ/T,−1/T

(62)

Properties of diffusion matrix

D= (D_ij), D_ij = Z

B

E^i+jv⊗d₀ dk

4π³ ∈R^3×3 and d₀ solves Q_el(d₀) =−F_µ,T(1−F_µ,T)v(k)

Proposition

D symmetric: D₀₁=D₁₀ and D_ij^>=D_ji

If(d₀,E(k)d₀)linearly independent then Dpositive definite Proof:

Symmetry: follows from symmetry of Q_el Show that forz ∈R⁶,z 6= 0,

z^>Dz = 1 2

Z

B²

σ_el(k,k⁰)δ(E⁰−E)

z· d0

Ed₀

2 dk⁰dk

4π³F(1−F) >0 since z·(d0,Ed0)^>= 0 would imply linear dependence of (d0,Ed0).

(63)

Properties of relaxation-time term

Inelastic (electron-phonon) collision operator:

Qin(f) = Z

B

s(k⁰,k)f⁰(1−f)−s(k,k⁰)f(1−f⁰) dk⁰ s(k,k⁰) =σ (1 +N)δ(E⁰−E+Eph) +Nδ(E⁰−E−Eph) where N: phonon occupation number,E_ph: phonon energy

Proposition

W is monotone, W(µ,T)(T −1)≤0 for all µ∈R, T >0 Proof: After some manipulations,

W(µ,T)(T −1) = Z

B²

(1−F)(1−F⁰)δ(E −E⁰+E_ph)E_phNe^{−(E−µ)/T}

× e^E^ph^/T−e^E^ph

(T −1)dk⁰dk 4π³ ≤0

(64)

Boundary conditions

∂tn−divJ0 = 0, ∂t(ne)−divJ1+∇V ·J0=W, x ∈Ω, t >0 J₀ =D₀₀

∇µ T −∇V

T

−D₀₁∇1

T, J₁=D₁₀

∇µ T −∇V

T

−D₁₁∇1 T Dirichlet conditions at contacts ΓD:

n=nD, T =TD, V =VD on ΓD

Neumann cond. at insulating boundary Γ_N: J₀·ν=J₁·ν =∇V ·ν = 0 on Γ_N

G

_D

G

_D

G

_N

W ^G

^D

G

_N

Improved boundary conditions for drift-diffusion (Yamnahakki 1995):

n+αJ0·ν =nD on ΓD

(second-order correction from Boltzmann equation)

Open problem: improved boundary conditions for energy-transport

(65)

Explicit models: spherical symmetric energy band

Assumptions:

Fµ,T approximated by MaxwellianM = exp(−(E−µ)/T) Scattering rate: σel(x,k,k⁰) =s(x,E(k)) for E(k) =E(k⁰)

Energy band spherically symm. monotone,|k|²=γ(E(|k|)), k ∈R³ Proposition

n ne

= e^µ/T 2π²

Z ∞ 0

e^−ε/Tp

γ(ε)γ⁰(ε) 1

ε

dε D_ij = e^µ/T

3π³ Z ∞

0

e^−ε/T γ(ε)ε^i+j

s(x, ε)γ⁰(ε)²dε, i,j = 0,1 Proof: Use coarea formula for, for instance,

n= Z

R³

e−(E(|k|)−µ)/T dk 4π³ = 1

4π³ Z ∞

0

Z

{E(ρ)=ε}

(...)dS_εdε

(66)

Parabolic band approximation

Assumptions:

Energy band: E(k) = ¹₂|k|²,k∈R³ Scattering rate: s(x, ε) =s1(x)ε^β,β≥0 Proposition

n =NT^3/2e^µ/T, N= 2

(2π)^3/2 density of states, ne= 3 2nT D=C(s1)Γ(2−β)nT^1/2−β

1 (2−β)T (2−β)T (3−β)(2−β)T²

Proof: Since γ(ε) = 2ε, n=

√2 π² e^µ/T

Z ∞ 0

e^−ε/T√ εdε=

√2

π²e^µ/TT^3/2Γ(³₂) = 2

(2π)^3/2T^3/2e^µ/T

(67)

Parabolic band approximation

Scattering rate: s(x, ε) =s1(x)ε^β,β ≥0 Diffusion matrix: typical choices for β

β = ¹₂ : Chen model D=

√π

2 C(s₁)n 1 ³₂T

3

2T ¹⁵₄T²

!

β = 0 : Lyumkis modelD=C(s₁)nT^1/2

1 2T 2T 6T²

Relaxation-time term:

W =−3 2

n(T −1)

τβ(T) , τ_β(T) =C(β,s₁)T^1/2−β Chen model: τβ constant inT

Kinetic and Macroscopic Models for Semiconductors