# A Quasi-Two-Dimensional Physics-Based Model of HEMTs without Smoothing Functions for Joining Linear and Saturation Regions of I-V Characteristics

Authors: Eugeny A. Ryndin; Amgad A. Al-Saman; Boris G. Konoplev
Journal: Active and Passive Electronic Components
 (2019)

---

## Abstract

A quasi-two-dimensional physics-based model of HEMT transistor without using any smoothing functions for joining the linear and saturation regions of current-voltage (I-V) characteristics was developed. Considering the intervalley transitions of electrons and the presence of holes in the channel of transistor, we calculated the nonuniform spatial distributions of the electrical field, electron temperature, and electron mobility within the channel. The model is in a good agreement with experimental data over the linear and saturation regions of operation. The model provides precise simulating of HEMT transistors and can be utilized as a tool for analysis and prediction of influence of the material parameters on device and circuit characteristics.


---

## Body

## 1. Introduction

High-electron-mobility transistor (HEMT) has appeared as a promising candidate for high-speed microwave and RF applications [1, 2]. High interest in the HEMT transistor and its modelling may be associated with its high-electron mobility, along with high-electron sheet concentration and excellent performance at very high frequency [3, 4]. To employ the full potential of these devices, accurate and fast simulation of the circuit based on HEMT is required, which needs compact models simulating their output characteristics. The output characteristics of HEMT transistors can be simulated using several modelling techniques: electrical (semiempirical, empirical) models or physics-based models. Despite the high accuracy of two- and three-dimensional physics-based numerical models with the rigorous field-dependent relation of carrier mobility within the channel, they are not convenient for utilization in the circuit design programs due to their high complexity [5–8].Recently reported semiempirical models show a good agreement with experimental data; nevertheless, they include many fitting parameters, which have nonphysical meaning [9–12]. This leads to remarkable limitations in their practical usage. With these difficulties, for example, complexity and time-consuming computations, high-speed optimized models based on device physics become preferable. A measure of proposed physics-based models of HEMTs has encountered severe limitations [13–15]. This is because they are usually developed using regional approximations that they are joined together by proper smoothing functions over the moderate region of the I-V characteristic of HEMT transistor.Attempts to overcome these considerable difficulties using surface potential-based models have been made by [16–20]. Although this approach is interesting, it fails to take into account two critical factors: (1) the dependence of carrier mobility on field spatial distribution and (2) the electronic temperature of 2-DEG within the channel. Moreover, to the best of the authors’ knowledge, the presence of holes in the channel was not incorporated into the charge conservation law by currently available physics-based models. Thus, this article aimed to introduce a novel quasi-two-dimensional physics-based approach for modelling HEMT transistors without using any smoothing functions for joining the linear and saturation regions of I-V characteristics. Furthermore, the field dependence of electron mobility will be considered as well as the intervalley transitions of electrons, the change of electron temperature, and the presence of holes within the channel.


## 2. Model Description

The proposed model has been developed for HEMT structure whose cross-section is illustrated in Figure1. It consists of a semi-insulating GaAs substrate on which an undopedi-GaAs layer of thickness H is grown, followed by a spacer layer of undopedi-Al1-xGaxAs of thickness Ws and a barrier layer of dopedn-Al1-xGaxAsof thickness WB with donor concentration ND. The highly doped ohmic contact regions in this structure are self-aligned with the gate.Figure 1
Structure of the AlGaAs/GaAs HEMT transistor.


### 2.1. Drain Current Calculation

Calculation of I-V characteristics of HEMT transistors is considered as one of the most significant objectives of modelling semiconductor devices. To simulate their I-V characteristics, the axesx andy would be pointed perpendicular and along the direction of current flow within channel correspondingly. Thus, drain current may be obtained by integration of electronic current density Jnx,y over the cross-section of the channel:(1)Id=∬SJnx,yⅆS=W∫0HJnx,yⅆx,where H is a thickness of the channel and W is the gate width.In our case, only the longitudinal drift component of electron current density will be taken into account. Consequently, (1) can be rewritten as(2)Id=WEyμyq∫0Hnx,yⅆx=-WEyμyQny,where E(y) is the distribution of electric field within the channel, μ(y) is the electron mobility, n(x,y) is the electron concentration, and Qn(y) is the surface charge density of electrons within the channel.In case of undoped channel HEMT transistor, besides the electrons depleted from barrier layer, there are intrinsic carriers (electrons and holes). Their densities are determined by the generation and recombination processes. It is necessary to notice that the drain current in such devices is only due to the electrons. Nevertheless, in place to calculate the surface charge density of electrons within the channel adequately, we should consider the presence of the positive charge within depletion region at the drain end. In case of MESFET, increasing of positive charges at depletion region is substantially attributed to ionization of immobile donors along with a nonsignificant increase in hole concentration, while in HEMT transistor this increase is associated only with hole concentration growing-up because of the undoped channel of HEMT. Namely, the growth in holes’ concentration should be reflected in the charge conservation law:(3)Qny+Qpy+QG+QB=0,where Qp(y) is the holes surface charge density, QG is the gate surface charge density, and QB is the surface charge density in the barrier layer.The last two variablesQG and QB in (3) can be expressed in terms of a gate voltage Ug, the barrier layer width (WB+WS), a donor concentration in the barrier layer ND, and the potential of channel Uy as follows:(4)QGUg,Uy=εε0WB+WSUg-Uy;(5)QBND,Ug=qNDWB,where ε is the relative permittivity of the barrier layer and ε0 is the permittivity of free space. The parameter’s values are given in Table 1 [22].Table 1
List of parameters used in simulation [22].
Symbol Physical meaning Value Unit ε The permittivity of the barrier layer 12.1 - ε 0 The permittivity of free space 8.85 × 10−12 F/m q Elementary charge 1.62 × 10−19 C k B The Boltzmann constant 1.38 × 10−23 J/K h The Planck constant 6.62 × 10−34 J⋅s T 0 The room temperature 300 K T L The lattice temperature 300 K E g , 0 The parameter of band gap model 1.52 eV α The parameter of band gap model 5.58 × 10−4 eV/K β The parameter of band gap model 220 K μ n , Γ m i n Electron mobility of highly doped GaAs at theΓ-valley 800 cm2/(V⋅s) μ n , X m i n Electron mobility of highly doped GaAs at theX-valley 40 cm2/(V⋅s) μ n , Γ , T 0 L Lattice mobility (T0 = 300 K) for GaAs at the Γ-valley 8500 cm2/(V⋅s) μ n , X , T 0 L Lattice mobility (T0 = 300 K) for GaAs at the X-valley 410 cm2/(V⋅s) C I Concentration of ionized impurity in the channel 1 × 1014 cm−3 С n r e f Temperature dependent fitting parameter 1 × 1017 cm−3 α n Temperature dependent fitting parameter 0.5 - γ 0 Fitting parameter -2.2 - β n Fitting parameter 2 - v s a t , n , Γ The electrons saturating velocity in theΓ-valley 2.5 × 107 cm/c v s a t , n , X The electrons saturating velocity in theX-valley 0.9 × 107 cm/c Δ The electronic gap betweenΓ and X valleys 0.31 eV τ w , 0 Energy relaxation model parameter for GaAs 0.48 ps τ w , 1 Energy relaxation model parameter for GaAs 0.025 ps C τ , 1 Energy relaxation model parameter for GaAs -0.053 - C τ , 2 Energy relaxation model parameter for GaAs 0.853 - C τ , 3 Energy relaxation model parameter for GaAs 0.5 - m n , Γ ∗ The effective mass of electrons inΓ-valley 5.92 × 10−32 kg m n , X ∗ The effective mass of electrons in X-valley 4.29 × 10−31 kg m p ∗ The effective mass of holes 4.19 × 10−31 kgSubstituting (4) and (5) into (3) yields(6)Qny+Qpy+εε0WB+WSUg-Uy+qNDWB=0.Considering the real HEMT transistors, their current-voltage characteristic is significantly affected by the charge of surface states. For that reason, the charge of the surface statesNs (its value varies from one device to another and must be extracted from the results of measurement) must be included into the charge conservation law:(7)Qny+Qpy+εε0WB+WSUg-Uy+qND+NsWB=0.With the assumption that the gate leakage current is negligible, a quasi-equilibrium should be established along the perpendicular direction to the channel plane so that we have(8a)np=ni2;(8b)ni=NcNv·exp⁡-Eg2kBTL;(8c)NC=22πmn,Γ∗kBTLh23/2;(8d)Nv=22πmp∗kBTLh23/2;(8e)Eg=Eg,0-α∙TL2β+TL,where ni is the intrinsic concentration, p is the hole concentration, Nc and Nv are the effective density of states in the conduction and valance bands, respectively, kB is the Boltzmann constant, TL is the lattice temperature, h is the Planck constant, T0 is the room temperature, Eg is the temperature-dependent band gap of channel semiconductor material Eg,0, α and β are the parameters of band gap model, mn,Γ∗ is the effective mass of electrons in Γ valley, and mp∗ is the effective mass of holes [22].The corresponding values of mentioned parameters used in simulation are listed in Table1 [22].The holesQp(y) and electron Qny surface charge densities can be expressed in terms of hole p(x,y) and electron n(x,y) concentrations using the following equations:(9)Qpy=q∫0Hnx,yⅆx≈qHp-y;(10)Qny=-q∫0Hnx,yⅆx≈-qHn-y,where p-(y) and n-(y) are hole and electron concentrations averaged over the channel thickness, respectively.Having assumed that replacement ofn and p in (8a) by p-(y) and n-(y) does not change the validity of (8a), it can be rewritten as(11)n-yp-y=ni2.Using expressions (9) and (10), (11) may be rewritten as follows: (12)QpyQny=-q2H2ni2.ExtractingQpy from (12) followed by substituting it into (7) leads to the following quadratic equation: (13)Qn2y+QnyqND+NsWB+εε0WB+WSUg-Uy-q2H2ni2=0.In order to determineQn(y), the roots of the quadratic equation (13) must be extracted so that we have (14)Qny=-εε02WBUsat-Uy±Usat-Uy2+b2.In (14), we have assumed that WB+WS≈WB and Usat and b can be calculated using the following expressions:(15)Usat=WBεε0qND+NsWB+Ug;b=2qHniWBεε0.By omitting the root of the quadratic equation with a negative sign, as it has no physical meaning, followed by substitutingQny into (2), one can derive (16)Id=Wεε02WBEyμyUsat-Uy+Usat-Uy2+b2.In (16), the spatial distributions Ey, Uy, and μy are interrelated functions of y, especially when one cannot ignore the facts that (1) the electron mobility depends on the field distribution and the electronic temperature within the channel and (2) the channel potential is distributed nonuniformly. To clarify these interrelationships, an adequate model of electron mobility and accurate calculation of nonuniformly distributed channel potential are required.


### 2.2. Electron Mobility

Theoretically, the electron mobility in semiconductors depends on the band structure of material concerned and on several scattering mechanisms. An accurate calculation of the band structure and distinct scattering mechanisms is extremely difficult and time-consuming. To circumvent this difficulty, the number of approximations was reported [23–25]. The electron mobility dependence on impurity concentration and lattice temperature may be expressed using formula proposed by Caughey and Thomas [22]:(17)μn,vLI=μn,vmin+μn,vL-μn,vmin1+CI/Сnrefαn,v=Γ,X,L,where CI is the concentration of ionized impurity, μn,vL is a lattice mobility, μn,vmin, Сnref, and αn are temperature-dependent fitting parameters, and v denotes considered valley (see Table 1).The lattice mobilityμn,vL in (17) is a function of the lattice temperature and can be modelled by power law as follows:(18)μn,vL=μn,v,T0LTLT0γ0,where TL is the lattice temperature, which in current study was considered as coordinate-independent and takes value equal to the room temperature, μn,v,T0L is the electron mobility at room temperature (i.e., at T0 = 300 K), and γ0 is a fitting parameter (see Table 1).With increasing electric field within the channel, it is necessary to allow for the impact of electrical field on the electron mobility. A measure of model options considering the reduction in the electron mobility caused by increasing electric field is available. The commonly used model is given by [22]:(19)μn,νLIFFny=μn,vLI1+μn,vLI∙Fny/vsat,n,vβn1/βn,where βn is a fitting parameter, vsat,n,v is the saturating velocity of electrons (see Table 1), and μn,vLI is the electronic mobility given by Caughey and Thomas.The driving forceFn in (19) can be written as(20)Fn=grad⁡Uy-1ngrad⁡kB·TLq·n.Assuming that the contribution of second term in the right-hand side of (20) is negligible compared to the first one, that is, 1/ngrad⁡kB·TL/q·n≪grad⁡U(y), we can rewrite (20) to take a form: (21)Fn≈grad⁡Uy.Equation (19) describes the field dependence of electron mobility; however, along with a decline in the electron mobility, increasing electric field within channel results in electron’s heating-up, especially near to the drain. When the electron temperature growth Tn(y) becomes comparable with the electronic gap between valleys Δ, the intervalley transitions should be taken into consideration [26].The decrease in mobility, in this case, depends on the ratio ofΓ- and X-valley populations, which can be found by [26]: (22)nXnΓ=4mn,X∗mn,Γ∗3/2exp⁡-ΔkBTny,where nX is the concentration of electrons in X-valley, nΓ is the concentration of electrons in Γ-valley,  kB is the Boltzmann constant, and mn,X∗ and mn,Γ∗ are the effective mass of electrons in X- and Γ-valleys, respectively (see Table 1).The concentration of electrons within the channel of HEMT transistor can be expressed as a sum of electron concentrations in the lower and upper valleys (n=nΓ+nX). Consequently, the overall electron mobility may be expressed as(23)μy=nΓμn,ΓLIF+nXμn,XLIFnΓ+nX,where μn,XLIF,μn,ΓLIF are the high-field electron mobility in X- and Γ-valleys correspondingly.Dividing the denominator and numerator on the right-hand side of (23) by nΓ, followed by applying (22), one can obtain the following expression for electron mobility with allowance for the intervalley transition:(24)μy=μn,XLIFRexp⁡-Δ/kBTny+μn,ΓLIF1+Rexp⁡-Δ/kBTny,where R is the state concentration ratio between X-valley and Γ-valley determined by the following equation: (25)R=4mn,X∗mn,Γ∗3/2.The spatial distribution of electronic temperatureTny in (22) and (24) can be calculated in terms of energy relaxation time τw, driving force Fn, and mobility μy as follows [22]: (26)Tny=TL+2q3kBτwμyFn2;(27)τw=τw,0+τw,1∙exp⁡Cτ,1∙TnyT02+Cτ,2TnyT0+Cτ,3TLT0,where τw,0, τw,1, Cτ,1, Cτ,2, and Cτ,3 are the parameters of the energy relaxation time model (see Table 1).


### 2.3. Channel Potential Calculation

In contrast to earlier existing traditional methods, we will consider the nonuniform distribution of potential within the channelUy. For that end, we start calculating the channel potential with the assumption that the gate and substrate currents are extremely negligible. This results in unaltered current flow through the channel, which means that the derivative of the drain current with respect to y must be equal to zero. Applying this in (16), one can derive a second-order differential equation with conditions U0=0 and ULg=Ud at the boundaries as follows:(28)ddyUsat-Uy+Usat-Uy2+b2μydUydy=0.Apparently, from (28), the electron mobility and the channel potential interrelate and depend on y in such intricate way, which must be considered when one solves (28). Under this circumstance, searching for an analytical solution will be tricky. For this reason, we have chosen to solve (28) numerically using finite difference method. Finite difference approximation of the second-order differential equation (28) for one-dimensional grid G=yj∣j=1,2,…,J has the following form:(29)Uj-Usat-Uj-Usat2+b2μjUj+1-Ujyj+1-yj--Uj-1-Usat-Uj-1-Usat2+b2μj-1Uj-Uj-1yj-yj-1=0,where j is the index of grid points.It is worth noting that the discretization scheme used was chosen to achieve fast convergence and accurate simulation result.The numerical procedure involves the following steps: first, we applied along they-axis 1D uniform mesh G with step equal to L/(J-1), where L is the channel length and J is mesh node which in our case was taken equal to 200. In order to calculate the spatial potential distribution by (28) and its discretized form (29), one should determine the initial values for the potential along with electron mobility, energy relaxation time, and carrier temperature. To do so, we arbitrarily assumed that the electrostatic potential distributed linearly throughout the channel, so the initial value of potential in each grid point may be calculated by(30)Uj0=U10+UJ0-U10yj-y1yJ-y1,where Uj(0) is the initial value of potential at the jth-node, U1(0) and UJ(0) are the potential values at the boundary, and y1 and yJ are the first and last grid points, respectively.As soon the initial distribution of the potential was determined, we used it correspondingly to calculate the initial values of electron mobility, energy relaxation time, and carrier temperature at each node of the grid. It is clearly seen from (22)–(27) that these quantities are interrelated which pose some calculation difficulties. For that reason, we decided to solve them self-consistently using nested loop: first we calculated the field-dependent electron mobility by (19) substituting (17), (18), and (21); then with the initial value of electron temperature equal to 300 K, we iteratively calculated the electron mobility, energy relaxation time, and carrier temperature using(31)μjk=μn,XLIFjkRexp⁡-Δ/kBTnjk+μn,ΓLIFjk1+Rexp⁡-Δ/kBTnjk;(32)τwjk=τw,0+τw,1∙exp⁡Cτ,1∙TnjkT02+Cτ,2TnjkT0+Cτ,3TLT0;(33)Tnjk+1=TL+2q3kBτwjkμjk∙Fnk2,where k is the iteration number.Initial values calculated by (30)–(33) have been utilized into (29) which in turns has been solved iteratively. In each iteration, the electron mobility, energy relaxation time, and electron temperature have been modified using recently calculated potential distribution. The result of this procedure was the spatial distributions of channel potential, electronic temperature, and electron mobility along the channel for each fixed gate and drain voltages.Last can be applied to (16) in order to extract the I-V characteristic. Foregoing calculation of the drain current was made for self-aligned HEMT structures; otherwise, one should include the influence of the resistance of the channel regions that do not locate beneath the gate of the transistor on the drain current. For that end, it is vital to invoke the following boundary conditions at the source and drain side of the channel:(34a)Uch0=U0;(34b)UchLg=Ud-IdRd+Rs,where Rd and Rs are the drain and source resistances, respectively.By applying these boundary conditions to (28) for each subsequent iterative step, one can calculate the drain current allowing for drain and source resistances.


## 2.1. Drain Current Calculation

Calculation of I-V characteristics of HEMT transistors is considered as one of the most significant objectives of modelling semiconductor devices. To simulate their I-V characteristics, the axesx andy would be pointed perpendicular and along the direction of current flow within channel correspondingly. Thus, drain current may be obtained by integration of electronic current density Jnx,y over the cross-section of the channel:(1)Id=∬SJnx,yⅆS=W∫0HJnx,yⅆx,where H is a thickness of the channel and W is the gate width.In our case, only the longitudinal drift component of electron current density will be taken into account. Consequently, (1) can be rewritten as(2)Id=WEyμyq∫0Hnx,yⅆx=-WEyμyQny,where E(y) is the distribution of electric field within the channel, μ(y) is the electron mobility, n(x,y) is the electron concentration, and Qn(y) is the surface charge density of electrons within the channel.In case of undoped channel HEMT transistor, besides the electrons depleted from barrier layer, there are intrinsic carriers (electrons and holes). Their densities are determined by the generation and recombination processes. It is necessary to notice that the drain current in such devices is only due to the electrons. Nevertheless, in place to calculate the surface charge density of electrons within the channel adequately, we should consider the presence of the positive charge within depletion region at the drain end. In case of MESFET, increasing of positive charges at depletion region is substantially attributed to ionization of immobile donors along with a nonsignificant increase in hole concentration, while in HEMT transistor this increase is associated only with hole concentration growing-up because of the undoped channel of HEMT. Namely, the growth in holes’ concentration should be reflected in the charge conservation law:(3)Qny+Qpy+QG+QB=0,where Qp(y) is the holes surface charge density, QG is the gate surface charge density, and QB is the surface charge density in the barrier layer.The last two variablesQG and QB in (3) can be expressed in terms of a gate voltage Ug, the barrier layer width (WB+WS), a donor concentration in the barrier layer ND, and the potential of channel Uy as follows:(4)QGUg,Uy=εε0WB+WSUg-Uy;(5)QBND,Ug=qNDWB,where ε is the relative permittivity of the barrier layer and ε0 is the permittivity of free space. The parameter’s values are given in Table 1 [22].Table 1
List of parameters used in simulation [22].
Symbol Physical meaning Value Unit ε The permittivity of the barrier layer 12.1 - ε 0 The permittivity of free space 8.85 × 10−12 F/m q Elementary charge 1.62 × 10−19 C k B The Boltzmann constant 1.38 × 10−23 J/K h The Planck constant 6.62 × 10−34 J⋅s T 0 The room temperature 300 K T L The lattice temperature 300 K E g , 0 The parameter of band gap model 1.52 eV α The parameter of band gap model 5.58 × 10−4 eV/K β The parameter of band gap model 220 K μ n , Γ m i n Electron mobility of highly doped GaAs at theΓ-valley 800 cm2/(V⋅s) μ n , X m i n Electron mobility of highly doped GaAs at theX-valley 40 cm2/(V⋅s) μ n , Γ , T 0 L Lattice mobility (T0 = 300 K) for GaAs at the Γ-valley 8500 cm2/(V⋅s) μ n , X , T 0 L Lattice mobility (T0 = 300 K) for GaAs at the X-valley 410 cm2/(V⋅s) C I Concentration of ionized impurity in the channel 1 × 1014 cm−3 С n r e f Temperature dependent fitting parameter 1 × 1017 cm−3 α n Temperature dependent fitting parameter 0.5 - γ 0 Fitting parameter -2.2 - β n Fitting parameter 2 - v s a t , n , Γ The electrons saturating velocity in theΓ-valley 2.5 × 107 cm/c v s a t , n , X The electrons saturating velocity in theX-valley 0.9 × 107 cm/c Δ The electronic gap betweenΓ and X valleys 0.31 eV τ w , 0 Energy relaxation model parameter for GaAs 0.48 ps τ w , 1 Energy relaxation model parameter for GaAs 0.025 ps C τ , 1 Energy relaxation model parameter for GaAs -0.053 - C τ , 2 Energy relaxation model parameter for GaAs 0.853 - C τ , 3 Energy relaxation model parameter for GaAs 0.5 - m n , Γ ∗ The effective mass of electrons inΓ-valley 5.92 × 10−32 kg m n , X ∗ The effective mass of electrons in X-valley 4.29 × 10−31 kg m p ∗ The effective mass of holes 4.19 × 10−31 kgSubstituting (4) and (5) into (3) yields(6)Qny+Qpy+εε0WB+WSUg-Uy+qNDWB=0.Considering the real HEMT transistors, their current-voltage characteristic is significantly affected by the charge of surface states. For that reason, the charge of the surface statesNs (its value varies from one device to another and must be extracted from the results of measurement) must be included into the charge conservation law:(7)Qny+Qpy+εε0WB+WSUg-Uy+qND+NsWB=0.With the assumption that the gate leakage current is negligible, a quasi-equilibrium should be established along the perpendicular direction to the channel plane so that we have(8a)np=ni2;(8b)ni=NcNv·exp⁡-Eg2kBTL;(8c)NC=22πmn,Γ∗kBTLh23/2;(8d)Nv=22πmp∗kBTLh23/2;(8e)Eg=Eg,0-α∙TL2β+TL,where ni is the intrinsic concentration, p is the hole concentration, Nc and Nv are the effective density of states in the conduction and valance bands, respectively, kB is the Boltzmann constant, TL is the lattice temperature, h is the Planck constant, T0 is the room temperature, Eg is the temperature-dependent band gap of channel semiconductor material Eg,0, α and β are the parameters of band gap model, mn,Γ∗ is the effective mass of electrons in Γ valley, and mp∗ is the effective mass of holes [22].The corresponding values of mentioned parameters used in simulation are listed in Table1 [22].The holesQp(y) and electron Qny surface charge densities can be expressed in terms of hole p(x,y) and electron n(x,y) concentrations using the following equations:(9)Qpy=q∫0Hnx,yⅆx≈qHp-y;(10)Qny=-q∫0Hnx,yⅆx≈-qHn-y,where p-(y) and n-(y) are hole and electron concentrations averaged over the channel thickness, respectively.Having assumed that replacement ofn and p in (8a) by p-(y) and n-(y) does not change the validity of (8a), it can be rewritten as(11)n-yp-y=ni2.Using expressions (9) and (10), (11) may be rewritten as follows: (12)QpyQny=-q2H2ni2.ExtractingQpy from (12) followed by substituting it into (7) leads to the following quadratic equation: (13)Qn2y+QnyqND+NsWB+εε0WB+WSUg-Uy-q2H2ni2=0.In order to determineQn(y), the roots of the quadratic equation (13) must be extracted so that we have (14)Qny=-εε02WBUsat-Uy±Usat-Uy2+b2.In (14), we have assumed that WB+WS≈WB and Usat and b can be calculated using the following expressions:(15)Usat=WBεε0qND+NsWB+Ug;b=2qHniWBεε0.By omitting the root of the quadratic equation with a negative sign, as it has no physical meaning, followed by substitutingQny into (2), one can derive (16)Id=Wεε02WBEyμyUsat-Uy+Usat-Uy2+b2.In (16), the spatial distributions Ey, Uy, and μy are interrelated functions of y, especially when one cannot ignore the facts that (1) the electron mobility depends on the field distribution and the electronic temperature within the channel and (2) the channel potential is distributed nonuniformly. To clarify these interrelationships, an adequate model of electron mobility and accurate calculation of nonuniformly distributed channel potential are required.


## 2.2. Electron Mobility

Theoretically, the electron mobility in semiconductors depends on the band structure of material concerned and on several scattering mechanisms. An accurate calculation of the band structure and distinct scattering mechanisms is extremely difficult and time-consuming. To circumvent this difficulty, the number of approximations was reported [23–25]. The electron mobility dependence on impurity concentration and lattice temperature may be expressed using formula proposed by Caughey and Thomas [22]:(17)μn,vLI=μn,vmin+μn,vL-μn,vmin1+CI/Сnrefαn,v=Γ,X,L,where CI is the concentration of ionized impurity, μn,vL is a lattice mobility, μn,vmin, Сnref, and αn are temperature-dependent fitting parameters, and v denotes considered valley (see Table 1).The lattice mobilityμn,vL in (17) is a function of the lattice temperature and can be modelled by power law as follows:(18)μn,vL=μn,v,T0LTLT0γ0,where TL is the lattice temperature, which in current study was considered as coordinate-independent and takes value equal to the room temperature, μn,v,T0L is the electron mobility at room temperature (i.e., at T0 = 300 K), and γ0 is a fitting parameter (see Table 1).With increasing electric field within the channel, it is necessary to allow for the impact of electrical field on the electron mobility. A measure of model options considering the reduction in the electron mobility caused by increasing electric field is available. The commonly used model is given by [22]:(19)μn,νLIFFny=μn,vLI1+μn,vLI∙Fny/vsat,n,vβn1/βn,where βn is a fitting parameter, vsat,n,v is the saturating velocity of electrons (see Table 1), and μn,vLI is the electronic mobility given by Caughey and Thomas.The driving forceFn in (19) can be written as(20)Fn=grad⁡Uy-1ngrad⁡kB·TLq·n.Assuming that the contribution of second term in the right-hand side of (20) is negligible compared to the first one, that is, 1/ngrad⁡kB·TL/q·n≪grad⁡U(y), we can rewrite (20) to take a form: (21)Fn≈grad⁡Uy.Equation (19) describes the field dependence of electron mobility; however, along with a decline in the electron mobility, increasing electric field within channel results in electron’s heating-up, especially near to the drain. When the electron temperature growth Tn(y) becomes comparable with the electronic gap between valleys Δ, the intervalley transitions should be taken into consideration [26].The decrease in mobility, in this case, depends on the ratio ofΓ- and X-valley populations, which can be found by [26]: (22)nXnΓ=4mn,X∗mn,Γ∗3/2exp⁡-ΔkBTny,where nX is the concentration of electrons in X-valley, nΓ is the concentration of electrons in Γ-valley,  kB is the Boltzmann constant, and mn,X∗ and mn,Γ∗ are the effective mass of electrons in X- and Γ-valleys, respectively (see Table 1).The concentration of electrons within the channel of HEMT transistor can be expressed as a sum of electron concentrations in the lower and upper valleys (n=nΓ+nX). Consequently, the overall electron mobility may be expressed as(23)μy=nΓμn,ΓLIF+nXμn,XLIFnΓ+nX,where μn,XLIF,μn,ΓLIF are the high-field electron mobility in X- and Γ-valleys correspondingly.Dividing the denominator and numerator on the right-hand side of (23) by nΓ, followed by applying (22), one can obtain the following expression for electron mobility with allowance for the intervalley transition:(24)μy=μn,XLIFRexp⁡-Δ/kBTny+μn,ΓLIF1+Rexp⁡-Δ/kBTny,where R is the state concentration ratio between X-valley and Γ-valley determined by the following equation: (25)R=4mn,X∗mn,Γ∗3/2.The spatial distribution of electronic temperatureTny in (22) and (24) can be calculated in terms of energy relaxation time τw, driving force Fn, and mobility μy as follows [22]: (26)Tny=TL+2q3kBτwμyFn2;(27)τw=τw,0+τw,1∙exp⁡Cτ,1∙TnyT02+Cτ,2TnyT0+Cτ,3TLT0,where τw,0, τw,1, Cτ,1, Cτ,2, and Cτ,3 are the parameters of the energy relaxation time model (see Table 1).


## 2.3. Channel Potential Calculation

In contrast to earlier existing traditional methods, we will consider the nonuniform distribution of potential within the channelUy. For that end, we start calculating the channel potential with the assumption that the gate and substrate currents are extremely negligible. This results in unaltered current flow through the channel, which means that the derivative of the drain current with respect to y must be equal to zero. Applying this in (16), one can derive a second-order differential equation with conditions U0=0 and ULg=Ud at the boundaries as follows:(28)ddyUsat-Uy+Usat-Uy2+b2μydUydy=0.Apparently, from (28), the electron mobility and the channel potential interrelate and depend on y in such intricate way, which must be considered when one solves (28). Under this circumstance, searching for an analytical solution will be tricky. For this reason, we have chosen to solve (28) numerically using finite difference method. Finite difference approximation of the second-order differential equation (28) for one-dimensional grid G=yj∣j=1,2,…,J has the following form:(29)Uj-Usat-Uj-Usat2+b2μjUj+1-Ujyj+1-yj--Uj-1-Usat-Uj-1-Usat2+b2μj-1Uj-Uj-1yj-yj-1=0,where j is the index of grid points.It is worth noting that the discretization scheme used was chosen to achieve fast convergence and accurate simulation result.The numerical procedure involves the following steps: first, we applied along they-axis 1D uniform mesh G with step equal to L/(J-1), where L is the channel length and J is mesh node which in our case was taken equal to 200. In order to calculate the spatial potential distribution by (28) and its discretized form (29), one should determine the initial values for the potential along with electron mobility, energy relaxation time, and carrier temperature. To do so, we arbitrarily assumed that the electrostatic potential distributed linearly throughout the channel, so the initial value of potential in each grid point may be calculated by(30)Uj0=U10+UJ0-U10yj-y1yJ-y1,where Uj(0) is the initial value of potential at the jth-node, U1(0) and UJ(0) are the potential values at the boundary, and y1 and yJ are the first and last grid points, respectively.As soon the initial distribution of the potential was determined, we used it correspondingly to calculate the initial values of electron mobility, energy relaxation time, and carrier temperature at each node of the grid. It is clearly seen from (22)–(27) that these quantities are interrelated which pose some calculation difficulties. For that reason, we decided to solve them self-consistently using nested loop: first we calculated the field-dependent electron mobility by (19) substituting (17), (18), and (21); then with the initial value of electron temperature equal to 300 K, we iteratively calculated the electron mobility, energy relaxation time, and carrier temperature using(31)μjk=μn,XLIFjkRexp⁡-Δ/kBTnjk+μn,ΓLIFjk1+Rexp⁡-Δ/kBTnjk;(32)τwjk=τw,0+τw,1∙exp⁡Cτ,1∙TnjkT02+Cτ,2TnjkT0+Cτ,3TLT0;(33)Tnjk+1=TL+2q3kBτwjkμjk∙Fnk2,where k is the iteration number.Initial values calculated by (30)–(33) have been utilized into (29) which in turns has been solved iteratively. In each iteration, the electron mobility, energy relaxation time, and electron temperature have been modified using recently calculated potential distribution. The result of this procedure was the spatial distributions of channel potential, electronic temperature, and electron mobility along the channel for each fixed gate and drain voltages.Last can be applied to (16) in order to extract the I-V characteristic. Foregoing calculation of the drain current was made for self-aligned HEMT structures; otherwise, one should include the influence of the resistance of the channel regions that do not locate beneath the gate of the transistor on the drain current. For that end, it is vital to invoke the following boundary conditions at the source and drain side of the channel:(34a)Uch0=U0;(34b)UchLg=Ud-IdRd+Rs,where Rd and Rs are the drain and source resistances, respectively.By applying these boundary conditions to (28) for each subsequent iterative step, one can calculate the drain current allowing for drain and source resistances.


## 3. Results and Discussion

Our model (16)–(34b) was applied to simulate the I-V characteristic of HEMT transistor depicted in Figure 1 with a gate length Lg = 500 nm, H = 500 nm, Ws = 3 nm, WB = 45 nm, ND = 3 × 1017 cm−3, Ns = -5 × 1015 cm−3, and W = 2 μm. The core component of our model is to precisely compute the potential distribution for each point within the channel. This can be utilized to determine the electrical field distribution and the electron mobility, which in turn are required for the drain current calculation.Figure2 shows the numerically calculated spatial distributions of the potential (a), the electron mobility (b), the electron temperature (c), the electrical field (d), the electron drift velocity (e), and X- and Γ-valley population ratio (f) within the channel of HEMT transistor for different values of the drain voltage.Figure 2
Numerically calculated spatial distributions of the potential (a), the electron mobility (b), the electron temperature (c), the electrical field (d), the electron drift velocity (e)و and X- andΓ-valley population ratio (f) within the channel with the different value of the drain voltage Ud. The gate length and gate voltage are Lg = 500 nm and Ug = 0.5 V, respectively.
(a) (b) (c) (d) (e) (f)By the results depicted in Figure2, it was found that, with the gate length Lg = 500 nm and the drain voltage up to Ud = 0.1 V, the distributions of the electrical field are linear within the channel. However, while the drain voltage becomes more than 0.1 V, the electrical field gets increased within the channel and shows more nonuniformity. That in turn leads to an increase in the electron temperature and correspondingly decrease in electron mobility along with the electron drift velocity. As illustrated in Figure 2, the electronic temperature reaches its maximum value at the drain end as well as the decrease in mobility of electrons, which associated with the resistance increasing of the region around the drain. The increase in resistance leads to an intensification of electrical field distribution nonuniformity as clearly seen from Figure 2.Current-voltage characteristics of AlGaAs/GaAs HEMT transistors with the gate lengthLg = 500 nm were calculated by (16) (see Figure 3) allowing for the obtained spatial distributions of the electrical field and electron mobility within the channel. Examination of Figure 3 highlights that the drain current saturates at the drain voltage Ud = 0.05 – 0.2 V.Figure 3
Current-voltage characteristics of AlGaAs/GaAs HEMT with 500 nm gate length.To assess the validity of the proposed model, a comparative study has been accomplished, comparing simulated I-V curves with experimental data from [21] for GaAs/AlGaAs HEMT transistor with the gate length Lg = 500 nm, the barrier layer of thickness WB = 49 nm, the gate width W = 200 μm, the doping concentration ND = 1 × 1017 cm−3, the charge of surface states Ns = -5 × 1015 cm−3, and an undoped channel layer of thickness H = 500 nm.The results presented in Figure4 show that the agreement between experimental data and simulation results is good, which confirms the validity of proposed model to describe HEMT structures. Disagreement with experimental curves for small gate voltage can be attributed to the voltage dependence value of the source and drain resistances.Figure 4
Experimental [21] (dotted lines with circle marks) and simulation (solid lines) current-voltage characteristics of HEMT transistor with Lg = 500 nm.Having analyzed AlGaAs/GaAs HEMT transistors and providing a comparison with real AlGaAs/GaAs HEMT transistor, we have estimated the validity of the proposed model and its ability to predict I-V characteristics taking into account the effect of intervalley transition. Incorporating the physical device geometry in the proposed model makes it more capable of simulating different devices beyond the HEMT transistors.


## 4. Conclusions

A new quasi-two-dimensional physics-based approach for simulating HEMT transistors has been developed. The important identity of the proposed approach is that the calculation of the I-V characteristics has been made by (16) without applying any smoothing function for joining the linear and saturation regions of I-V curves. The impact of intervalley transitions on electron mobility has been assessed. The model has been verified by comparison with experimental I-V curves of GaAs/AlGaAs HEMT transistor with the gate length Lg = 500 nm. A good agreement with experimental data was achieved.



---

*Source: 5135637-2019-04-01.xml*