Two-dimensional Modelling of HEMTs Using Multigrids with Quantum Correction

The two-dimensional multi-layered HEMT is modelled isothermally by solving the Poisson and current continuity equations consistently with the Schrodinger equation. A multigrid method is used on the Poisson and current continuity equations while the electron density is calculated at each level by solving the Schrodinger equation in onedimensional slices perpendicular to the layer structure. A correction factor is introduced which enables relatively accurate solutions to be obtained using a low number of eigensolutions. A novel method for discretising the current density which can be generalised to the non-isothermal case is described. Results are illustrated using a two layer AlGaAs-GaAs HEMT.


INTRODUCTION
This paper addresses the problem of solving the modelling equations for High Electron Mobility Transistors (HEMTs) in which the Schr6dinger equation is solved explicitly.A basic 4-layer structure consisting of layers of AluGal_uAs with different proportions of mole fraction u of alumi- nium content is shown in Figure with the contacts arranged along the x-axis and with the y-axis perpendicular to the layer boundaries.Such devices have been shown to operate up to 213 GHz [1].
Fast operation is achieved because the electrons near the layer interfaces (the "2-dimensional" electrons) are confined in narrow potential wells and suffer little scattering.The equations governing the modelling of the device are the Schr6din- ger, Poisson, current continuity, and the energy transport equations.In this paper we will apply some new techniques to the isothermal case, although they may be easily generalised to the non-isothermal case.Here, q is the magnitude of the electron charge, n is the total electron density, and e0 is the permittivity of free space.Writing EF=--qc) as the quasi Fermi level, q5 is fixed at the bias values on the contacts.
The values of b on the source and drain are given by putting the electron density equal to the doping density using the Fermi function expressions described below.We put b Eh-bb-O5 on the gate where qSb =-0.7 V is the built-in potential.
(ii) In the steady state we take V J O, J -q#nVd? (1.2) for the current density J. Generally the mobility may depend on the field.(iii) The Schr6dinger equation equations which model the devices [1--13].A special function (the C-function) has been devel- oped [2] which allows both the current continuity and energy transport equations to be discretised in a consistent manner.A multigrid method [3] has been developed and applied to a one dimensional model, and a more efficient method of including the eigensolutions of the Schr6dinger equation has been developed [4].
Figure 2 shows the band offset structure between different layers.This structure will depend on the proportion u of A1 in the layer AluGal_As.The equations to be solved are: for the solution of the electrostatic potential b.
--V" Vi -[-Vxc -}-Eh q)i Ai(. where b is the electrostatic potential and (i and ,i (i=0, 1,2,...) are the energy eigenfunctions (normalised) and eigenvalues respectively.The forms of the exchange correlation energy Vx,, effective mass m, relative permittivity er, conduction band discontinuity Eh, N + and T-are taken from Adachi [5].The form of the kinetic energy operator (-h2/Z)V (1/m)V is taken for its he- rmitian qualities [14, 15].The full two-dimen- sional solution of the Schr6dinger equation is very time consuming, and in practice is solved in one dimensional slices perpendicular to the layer structure (that is, in the y-direction).We therefore solve the one dimensional equation by imposing the boundary conditions .i(0)i (Y) 0 (i 0, 1,2,...).
The electron density n inside the quantum wells will be given in terms of the eigensolutions of this equation.The total electron density is given by n n2+ n3 where n2 (for the "2-dimensional" electrons) is the contribution from the sub-bands given by the solution of the Schr6dinger equation and n3 (for the "3-dimensional" electrons) is the bulk electron density.We first choose a maximum number L of eigensolutions of the Schr6dinger equation with which to work.Outside the potential well defined by ,"L--1 < Ec we have n2--0 and where Nc 2(2mkT/h2) .5.
Inside the well we must avoid the double counting of contributions, and thus we take (for the one dimensional Schr6dinger solution only) L-1 n2 Nc2 [i(y)[2 ln(1 + e('/kr)(eF-Ai)) (1.6)   i=0 where Nc2 4wmkT/h 2, and n3- (1.7) Inside the well the integral in Equation (1.7) involves two parameters (1/kT)(L_I-E) and  (1/kT)(EF-A_), and its evaluation requires considerable computing time.An approximation to this integral has been developed [4] which speeds up the calculation by allowing a relatively small number of eigensolutions to be used.Outside the well the electron density is given in terms of the Fermi function F(/2) only, and the approximation of Bednarczyk and Bednarczyk [16] is used for the evaluation.
In the following section we briefly describe this correction process and introduce a new method of writing the current density so that the Bernoulli function method may be applied to the solution of the current-continuity equation.The final section describes the results of a two-dimensional simula- tion.
(2.3) Hence the value n3 inside the well may be rapidly calculated using Equation (2.2) with the value of c given by Equation (2.3).It was also shown that if this correction factor is used (rather than taking its value to be zero) then the number of eigensolu- tions needed for the calculation of n2 may be reduced to give good accuracy.We now describe a novel discretisation of the current-continuity equa- tion.In the isothermal case we may write J -q#nVc/:, (2.4) where EF =--q, and the electron density is a function of , , x and y: n n(, , x, y).Firstly consider the x-component of Equation (2.4).
Keeping y constant we have O---O--n + -0-n + nx when it is clear which variables are being kept constant.Substitution of the expression for O/Ox into Equation (2.4) gives the first component of J as J(x) q# +n n On n Ox (2.5)In the interval (xi <_ x <_ xi+l,y yj) define the quantity f(x, yj) q-i+1/2,j In n Then in this interval we can approximate Equation (2.5) by Of On J(x) -x n + fl 0-'-- where a =_ q#n/n4 and /3-q#n/n 4. We invoke the Bernoulli function approach [2] to this discretisation by assuming that/3, a/3 and Of/Ox are constants in the interval.Then J(x)i+1/2,j hi (fi+l,j fi,j)/ ni+l,j-n(() (j+,,j f',j))ni, j i+1/2,j (2.6) where B(x) =_ x/(e x-1) and hi Xi + 1--Xi.The y- component may be similarly discretised by intro- ducing a function g(xi, y) which is defined in the interval (x xi, yj <y _< Yj+ 1).Assumptions similar to the above have been derived elsewhere [2] and as a brief indication of their validity, note that in the non-degenerate case we have a//3 =-q/kT when n is given in purely exponential terms.The above formulation is general in that the exact functional form of the electron density is not specified.Once the specification is given then the derivatives n, ne, nx and ny may be calculated.In the non-isothermal case, Equation (2.4) will not apply and the generalisation of Equation (2.5) will contain extra terms On/OT.The solution will then be given using the C-function approach [2] which generalises the Bernoulli function approach.
In order to reduce the computational effort in reaching a solution, a multigrid approach has been developed [3] in which the Poisson and current- continuity equations are solved using multigrids while the Schr6dinger equation is solved by a non- multigrid scheme at each level to provide the expression for the electron density n.This has increased the speed of the solution by almost five times in some cases.

3.. SIMULATION RESULTS
The model problem is that of a two-layer A10.3 Ga0.7As-GaAs HEMT.The ends of the source, gate and drain were at x 0.0 lam, 0.1 lam, 0.4 lam, 0.7 lam, 1.3 gm and 1.4 gm, the total depth of the device was Y=0.3 lam with the layer interface at y 0.1 gm.The doping of the A1GaAs and GaAs was taken as 5x 1023 m -3 and x 1019 m -3 respec- tively.High doping regions of 1.2x 10 24 m-3 were taken around the source an drain to simulate ohmic contacts, and a depletion value of 8x 1016 m -3 was taken on the contact edges between source and gate and gate and drain to pin Ec to a value at which the quantum well could form.The electron temperature was taken as a constant 300 K.The mobilities of the source (a) drain \0 0.1 FIGURE 3 Results of the simulation with Vgs -0.1 V and Vas 2.0 V, with all distances in gm, showing (a) the conduction band E., (b) the electron density n, and c) the second energy eigenfunction 1. A1GaAs and GaAs were taken as 0.2 m2V IS--1 and 0.8 m2V-ls -1 respectively.

QUANTUM CORRECTIONS IN HEMTs
Results shown in Figure 3 were derived for the case Vgs =-0.1 V, Vds 2.0 V.A multigrid method was used to solve the Poisson and current continuity equations.The Schr6dinger equation was solved on each grid by a non-multigrid method to provide the electron density.The eigenvalues were found using a QL algorithm with implicit shifts, ordered, and the first L eigenfunc- tions found by back substitution.The multigrid method described in [3] was used.Four uniform grids were taken with the coarsest grid being 40x25 -this was found to be the minimum size for rapid calculation of the exact solution on the coarsest grid without losing the quantum well structure in the y-direction.The coarse grid solution was found using Newton iteration with automatic stopping based on the sizes of the residuals of the Poisson and current continuity equations.It was found that an initial number of 20 iterations had to be performed on this solution without quantum correction to allow the quantum well structure to become established.Ten pre-and post-smoothing iterations were made on the finer grids.As found earlier [4] it was found that using only three eigenvalues (L 2) with the correction factor of Equation (2.3) provided comparable results with those of the benchmark case L=9 with no correction used, thus speeding the numerical solution.Further, solutions could be found for fixed values of Vgs and Vds without having to approach them in small increments. FIGURE