Simplified energy-balance model for pragmatic multi-dimensional device simulation

Solid-State Electronics Vol. 41, No. 11, pp. 1795-1802.1997 © 1997ElsevierScienceLtd. All rights reserved Printed in Great Britain

Pergamon

P I h PII: S0038-1101(97)00142-1

0038-1101/97 $17.OO+ 0.OO

SIMPLIFIED ENERGY-BALANCE MODEL FOR PRAGMATIC MULTI-DIMENSIONAL DEVICE SIMULATION D U C K H Y U N CHANG and JERRY G. FOSSUM Department of Electricaland Computer Engineering,University of Florida, Gainesville,FL 32611, U.S.A. (Received 1 November 1996; in revised form 29 March 1997)

To pragmatically account for non-local carrier heating and hot-carrier effects such as velocity overshoot and impact ionization in multi-dimensional numerical device simulation, a new simplified energy-balance (SEB) model is developed and implemented in FLOODS[16] as a pragmatic option. In the SEB model, the energy-relaxation length is estimated from a pre-process dri•diffusion simulation using the carrier-velocitydistribution predicted throughout the devicedomain, and is used without change in a subsequent simpler hydrodynamic (SHD) simulation. The new SEB model was verified by comparison of two-dimensional SHD and full HD DC simulations of a submicron MOSFET. The SHD simulations yield detailed distributions of carrier temperature, carrier velocity,and impact-ionization rate, which agree well with the full HD simulation results obtained with FLOODS. The most noteworthy feature of the new SEB/SHD model is its computational efficiency,which results from reduced Newton iteration counts caused by the enhanced linearity. Relative to full HD, SHD simulation times can be shorter by as much as an order of magnitude since larger voltage steps for DC sweeps and larger time steps for transient simulations can be used. The improved computational efficiencycan enable pragmatic three-dimensional SHD device simulation as well, for which the SEB implementation would be straightforward as it is in FLOODS or any robust HD simulator. © 1997 Elsevier Science Ltd Abstract

l. INTRODUCTION In contemporary deep-submicron devices, nonstationary, or non-local carrier transport phenomena make prominent impacts on device performance and reliability. Since they are fundamentally described by the carrier energy distribution, not by the local electric field, the conventional drift-diffusion (DD) model is not adequate for device design in spite of its maturity and computational efficiency. Monte Carlo (MC) analysis of microscopic carrier transport is considered as the most reliable simulation method since it solves the Bolzmann transport equation (BTE) directly and thus provides the detailed energy distribution[l]. However its statistical nature leads to extremely long computation time which limits its utility for engineering applications[2]. The hydrodynamic (HD) model, in which carrier dynamics is described by averages of energy and velocity, was investigated widely as an alternative to the MC method because of its capability of describing non-local carrier transport in less computing time[3-6]. However the overall accuracy of the HD model is questionable because of equivocation in the implicit physical modeling, e.g., in the characterization of momentum- and energy-relaxation times and their dependences on carrier energy[7], and the difficulties in properly incorporating non-parabolic energy bands[8], non-Maxwellian carrier distributions[9], and higher-order moments of the BTE

(e.g., heat flux)[10]. Also, it suffers from much longer simulation time relative to that of the DD model. Indeed users of contemporary HD numerical device simulators can wait hours for a simulation to be executed, and still have little confidence in the results. To obtain an adequate description of the average carrier energy more efficiently than via the full HD model, simplification of the energy-balance equation was proposed recently for post-processing use. The post-processing technique, in which carrier energy is calculated from the electric field distribution obtained from a DD simulation, was primarily used to predict the energy-dependent impact-ionization rate[11-14]. Slotboom et al.[12] simplified the energy-balance equation by assuming a constant (saturated drift) velocity throughout the device, and then used it to model the impact ionization along a single current path in MOSFETs. However the model does not provide detailed information about the two-dimensional distribution of the carrier energy within the devices. Recently, Agostinelli and Tasch[14] demonstrated an energy-dependent two-dimensional impact-ionization, or substrate current model for the simulation of submicron MOSFETs, in which the simplified one-dimensional energy-balance equation proposed by Cook[15] is applied to many current contours in order to generate a two-dimensional representation of carrier energy. However the average energy predicted by this post-processing technique is not correct in highly scaled devices since it decouples

1795

1796

Duckhyun Chang and Jerry G. Fossum

the electric field from the energy. For example in deep-submicron MOSFETs, velocity overshoot influences the electric field distribution for a given bias condition, and effectively defines a higher drain saturation voltage, which in turn defines a higher current and affects gate charge, or capacitance, as well. This effect further reduces the impact-ionization current because of the lower electric field near the drain. The significant electric field-energy coupling hence tends to invalidate the post-processing method of characterizing the substrate current in scaled MOSFETs. Furthermore, issues related to the implicit simplified modeling of the carrier velocity[12,15] in the energy-balance equation prohibit use of these methods in general device simulation. The simplifications are not appropriate for regions in which the net flow of the carriers is small. According to these models, the carriers are heated, even at equilibrium, because of large built-in electric fields, for example at the base-emitter junction of a scaled BJT. This is contrary to the definition of thermal (or quasi-) equilibrium. In such cases, the carriers are not heated because the drift current is balanced out by the diffusion current. For a MOSFET, this error can be more serious for all bias conditions. In the vicinity of the inversion layer, the noted simplifications predict high carrier energy caused by the large transverse electric field and field gradient; but there is no carrier heating resulting from that field. Hence, the inadequate modeling of the carrier velocity in the energy-balance equation prohibits reliable device simulation, and indeed could cause numerical errors in simulation as well. In this work, we develop a reliable simplified energy-balance (SEB) model which accounts for the non-local effects with reasonable accuracy. The new SEB model is formulated for practical multidimensional device simulation and is implemented in the two-dimensional tool FLOODS (Florida Objected-Orient Device Simulator)[16] as a pragmatic option. The physical models underlying the SEB formalism are also discussed. Non-local effects in the submicron bulk-silicon MOSFET, such as velocity overshoot and impact ionization, are analyzed based on simulation results, and are compared with those of the full HD model in FLOODS. The superior computational efficiency of the new model, relative to the HD model, which renders it pragmatic and attractive for three-dimensional and transient device simulation, is addressed as well.

2. SEB MODEL A N D I T S I M P L E M E N T A T I O N IN F L O O D S

The conventional HD model comprises the first three moment equations of the BTE with Poisson's equation. For steady-state conditions, the moment

equations for electrons are[17]; V.(nv.)+(Ro

v~=-#,

E+

- Go) = 0

q ~--+V

(l)

(2)

v.. Vw0 = - q(v,. E) - 1 V. (nkB Tnv, + Q) -- Wo - w o (3) n

which describe conservation of the carrier density, average momentum, and average energy, respectively. In (1)-(3), n is the electron concentration, ( R o - G , ) is the net recombination rate, v. is the average electron velocity, #. is the electron mobility, E is the electric field, T, is the average electron temperature, w, is the average electron energy, Q is the heat flux, and rw is the energy-relaxation time. In 3 I * 2 the energy-balance eqn (3), w, = ~kBT, + ~m v, and wo = 3 k s T t where TL is the lattice temperature. In the HD model, the system of eqns (1)-(3) for electrons and analogous ones for holes is solved simultaneously with Poisson's equation. Inclusion of the velocity v., defined in eqn (2) as - J , / q n where J. is the electron current density, in the energy-balance eqn (3) greatly increases the nonlinearity of the system, and thus the computation time of simulation. The post-processing models described in Section 1 ameliorate this problem, but as noted they exaggerate the electron energy and hence, for example, overly predict the impact-ionization rates, and they are not suitable for general device simulation. We propose a new simplified energy-balance (SEB) model with reliable accuracy, still maintaining the computational efficiency like the post-processing model proposed by Slotboom[12]. We first assume that the drift energy is only a small part of the total kinetic energy; hence the average electron energy is approximated as w. ~ 3kB7",. Under most conditions in MOS devices, this assumption is generally valid, though near the source region, for example in a 0.2 pm device, the drift energy can be nearly 15% of the total kinetic energy[l]. We next assume that the heat-flux term Q is negligible in determining the electron temperature. This term in eqn (3) represents the third moment of the BTE and is thus an additional unknown[10]. It is generally approximated by the Fourier law, but in large temperature-gradient regions, this approximation greatly overestimates the heat conduction and thus leads to non-physical results[17]. To resolve this problem, a more sophisticated model for the heat flux was presented[l]; however, numerical difficulties caused by it provide an additional computational burden. Because of this modeling complexity and the fact that Q is often relatively small, the heat-flux term is frequently and justifiably neglected in the energy-balance model[7, ! 2]. Now focusing on eqn (3), we assume that V.(nv.)~0, which, as eqn (1) shows, is equivalent to neglecting carrier recombination/generation in

1797

Simplified energy-balance model for pragmatic multi-dimensional device simulation regions where carriers are heated. Then eqn (3) simplifies to

[ Set Bias [

v .VT,+3 (Tn -2 q v,'E. z-~TL) - - 5--~B

DD Simulation

(4)

Further assuming that the effect of diffusion (or Vn) is negligible since the electrons are heated mainly by drifting, and expressing the electric field on right-hand side by the gradient of the conduction band edge Ec, which renders eqn (4) applicable for homostructures as well as heterostructures, we write our steady-state SEB equation as;

5

z.~

Eo + ro '~wn,~j+ X,(ro -- to), = 0

(6)

where l,j is the length of the box side bisecting the triangle edge between mesh points i and j, Ai is the box area around point i, and 2wn,~j is the electron energy-relaxation length defined on triangle edge (i-j) of the grid, which is perpendicular to l~i. The energy-relaxation length, first defined by Slotboom[12] with constant v, = v.... is defined here as ,~wn = ~Zw.Vn 5

(7)

where v. is the velocity along each triangle edge and is thereby ~patially varying. In the first term of e qTq~_~_Eo n and T. are expressed as (Ec.i + Ecj)/2 and ~/T..~T.j, respectively, on each triangle edge[16]. There are two basic numerical approaches to solve the discretized semiconductor transport equations. One is Gummel's method in which the equations are decoupled and solved sequentially, and the other is

a-o

j=t~

Fig. 1. A two-dimensional grid structure with triangular meshes. SSE 41/11-E

I

(v,,~,p)

I

Calculate ~.,,~. ~,~p

on each triangle edge

-t

SHD Simulation:

For numerical device simulation, eqn (5) can be spatially discretized using the box method. For the two-dimensional case on a triangular mesh structure as illustrated in Fig. 1, the discretization gives ~ol, j

-3

[

Poisson's Eq. (~/) Carrier Continuity Eqs. (n, p ) SEB Eqs. (Tn, Tp )

@

I Set New Bias I Fig. 2. Flowcharted SEB/SHD DC model algorithm for a specific bias.

Newton's method in which the equations are solved simultaneously, thus accounting for all of the coupling between state variables[18]. The Newton algorithm is very stable and the solution time is nearly independent of the bias condition, even if high-level injection obtains. For the HD model, it is very difficult to directly implement the full Newton method because of the large Jacobian size and the nonlinearity introduced by the energy-balance equation. To alleviate this difficulty, FLOODS uses the Newton approach, but with the drift-diffusion simulation as a pre-processor. From the DD simulation, which uses Tn from the previous bias-point solution, the solver defines the initial guesses for variables of the subsequent HD simulation. This approach was verified as a reliable method to characterize T, > TL and the performance of scaled devices[16], albeit less computationally efficient than the decoupled Gummel method. Our new SEB model was developed based on the numerical HD method used in FLOODS. The key assumption is that the coupling between the carrier energy and the electric field is more significant than the one between the carrier energy and carrier velocity, which is reasonable because the carrier velocity tends to saturate with the carrier energy. The algorithm of the DC model is flowcharted in Fig. 2. The coupling between the electric field and carrier energies, which is strong in scaled devices, is faithfully accounted for, in contrast to post-processing

1798

Duckhyun Chang and Jerry G. Fossum

analytical methods for estimating carrier tempera- given by tures, e.g., the one proposed by Slotboom[12]. According to Slotboom's model, the electron G, = ~ [ ~ / ~ e x p ( - l ) - e r f c ( ~ l "~l (8, energy-relaxation length is assumed to be spatially constant, empirically set at 65 nm. However, we found that even though this assumption has some where @ - ks T./Wth and ro is a characteristic time merit for predicting the peak impact-ionization rate which depends on the effective mass of the in simple one-dimensional devices, it is not physical conduction-band electrons and the dielectric constant and it is not sufficiently accurate for use in general of the semiconductor. In eqn (8), the threshold e n e r g y w,h and the time constant z0 are treated as adjustable multi-dimensional device simulation. Alternatively in our SEB model, the energy-relax- parameters tuned via measurement. In FLOODS, ation length in eqn (7) is estimated from a pre-process this carrier generation caused by the impact drift-diffusion simulation using the carrier velocity ionization is implemented as a part of the distribution predicted throughout the device domain, carrier-continuity eqn (1). Alternatively, for computational efficiency, it can be calculated in post as can be seen in Fig. 2. In the preliminary DD simulation, the complete velocity expression (2) is processing, which is adequate for the common case of used, with the VTo term evaluated from the previous weak impact ionization, in contrast to the case of bias-point solution. The estimated values of 2w are avalanche breakdown, which necessitates the former then used without change in a subsequent simpler exact calculation. The substrate current of a hydrodynamic (SHD) simulation where Poisson's MOSFET is calculated by the integration of the equation, the carrier continuity equations with the generation rate over the simulation mesh. We stress velocity expressions, and the carrier SEB's, given here that our SEB/SHD model enables simulation of in eqn (5) for electrons, are solved simultaneously. impact ionization as well as other hot-carrier effects The carrier energies are thus directly coupled to directly in the coupled numerical iteration process. and defined by the electric field distribution; the The SEB/SHD model is not merely decoupled post pre-defined energy-relaxation length distribution processing. in the SEB equation is representative, yet eliminates dependences on the carrier densities, thereby 3. MODEL VERIFICATION substantially enhancing the linearity of the system The SEB/SHD model option was used in d.c. of SHD equations. The SEB/SHD model is reasonably accurate, and it can substantively increase simulations of the 0.5 #m bulk-silicon n-channel MOSFET used by Liang[16] for FLOODS demoncomputational efficiency because of the enhanced stration. The gate oxide thickness is 15 nm, the linearity. The SEB/SHD DC model was implemented as a source/drain junction depth is 0.3pro, and the pragmatic option in FLOODS without substantial channel doping density is 1.4 x 1017cm -3. An LDD code modification because of the high portability of the object-oriented device simulator[16]. The SHD algorithm and the discretization scheme used 6000 are the same as those implemented in the HD model H of FLOODS, but the original energy-balance 'O'--O SHD 5000 equation is replaced by the SEB equation, i.e., eqn (6) with 2,n in eqn (7) evaluated from the preceding DD simulation. The Jacobian elements are 4000 calculated accordingly. Note that the Jacobian for SHD simulation is much simpler than that of the conventional HD simulation due to the loosely 300O coupled nature of the carrier energies and densities. All physical models such as energy bands, SRH and Auger recombinations, impact-ionization rates, 2000 and energy-relaxation times were used without change. The energy dependence of the carrier 1000 mobility is incorporated through Klaassen's model[19]. Many expressions for the impact-ionization rate exist in the literature[ll,13,20]. In FLOODS, the 0.0 o.2 o.4 0.8 o.8 1.0 1.2 Scholl-Quade model[21] was used to describe it. The Dist~ee ( ~ ) model is derived from the BTE under the assumption that the spherical part of the electron distribution Fig. 3. Simulated electron temperature distribution along the 0.5 #m MOSFET channel for Vcs = VDs= 5.0 V. The function can be represented by a Maxwellian form. SHD and HD predictions of the location and magnitude of The resulting net generation rate per unit volume is the peak electron temperature show close agreement.

\,/77J

Simplified energy-balance model for pragmatic multi-dimensional device simulation 2.0

o-----o 81-ID ,-.

1.5 .

.,~

1.0

0.5

0.0 0.0

Ik~

0.2

0.4

0.6

0.8

1.0

.

1.2

Distance ( lan ) Fig. 4. Predicted average electron velocity distribution along the MOSFET channel corresponding to the simulations of Fig. 3. The SHD and HD models show consistent velocity overshoot near the drain, in contrast to the DD model. with 0.1 #m depth and 1.0 x 10~Scm-3 doping density is part of the device structure. Figure 3 exemplifies the electron temperature distribution along the channel simulated with the SEB/SHD model option in FLOODS and compares it with the full HD simulation result. To facilitate the direct comparison between the two models, equal values were used for model parameters like energy-relaxation time (Zw = 0.35 ps). Note that the strong non-equilibrium behaviours are predicted near the drain end of the channel by both models. The hoi-carrier temperature distribution simulated with the SEB/SHD model shows good agreement with the one from the HD model. Figure 4 shows the corresponding electron velocity obtained from the SHD and full HD simulations; the velocity, evaluated as vn = - J n / q n , was averaged over the channel depth. The SHD and HD velocity predictions are in close agreement. The velocity distribution obtained from the DD simulation is also included in the figure to emphasize the overshoot that actually occurs in the 0.5 gm device. The velocity attained by the electrons is seen to exceed its assumed (in FLOODS) s a t u r a t e d drift value of 9.8 x 10~cm s -j near the drain. This overshoot is produced primarily by an increase in electron mobility caused by Tn lagging E along predominant flowlines near the drain, and it becomes more significant as the channel length decreases. Note that the velocity overshoot can be beneficial to the performance of the device since it increases the drain current and reduces the sourcedrain transit time, but perhaps more significantly tends to reduce Tn for a given bias condition as we discussed in Section 1.

1799

Figure 5 shows the two-dimensional impact-ionization generation rate simulated with the SEB/SHD and full HD models at V~s = VDs = 5.0 V. The values used for the parameters W,h and r£ ~ in (8) are 4.0 eV and 7 x 10 t5 s -~, respectively. As seen in eqn (8), the impact-ionization rate is an increasing function of the carrier energy, and hence its distribution follows the temperature distribution. The SHD and HD predictions agree quite well, thereby verifying the validity of the SHD option. These two-dimensional generation rates were calculated via the post-processing option, which is adequate for the common case of weak impact ionization. As we mentioned in the previous section however, eqn (8), inserted in eqn (1), could be an integral part of the system. The substrate current of the MOSFET is calculated by the integration of the generation rate over the simulation mesh. For SHD and HD simulations, the substrate currents are 3.4 x 10 - s A p m -~ and 3.0 x 10-SA/~m -~, respectively, further confirming the validity of the SHD option. To demonstrate the real utility of the SHD option, we discuss the computational efficiency of the new model. Since the SHD simulation solves five nonlinear partial differential equations, i.e., Poisson's equation, the electron and hole continuity equations (with the carrier-velocity expressions), and the energy-balance equations (SEBs) for electrons and holes, its Jacobian has the same order as that of the full HD model. But, as we described in Section 2, the SHD model assumes fixed velocity, or energy-relaxation length distributions in the evaluation of the SEB equations. This simplification removes the n and p dependences of carrier energies and simplifies the Jacobian accordingly; the linearity of the system is thus enhanced. This enhanced linearity of the SEB/SHD model relative to the HD model results in reduced numerical iteration counts for a given bias, and hence gives shorter run times. Figure 6 shows the average iteration counts for the SHD and full HD DC MOSFET simulations done by sweeping Vos from 0 to 5.0 V in equal voltage steps A VDs, for different A Vos. As can be seen in the figure, for A VDs = 0.2 V, the SHD and HD simulations need 100 ( = 2 5 x 4) and 125 (=25 x 5) iterations, respectively, to obtain electron temperature distributions at VDs = 5.0 V; while for AVos --- 0.5 V, the SHD and HD simulations need 50 ( = 10 x 5) and 150 ( = 10 x 15) iterations, respectively. Note that as the simulation voltage step A VDs increases, the relative computational efficiency of the SEB/SHD model increases; that is, for the small voltage step, the simulation time of the SHD model is just 20% shorter than that of the full HD model, whereas for the large voltage step, the SHD simulation time is reduced by a factor of three relative to the HD simulation. This improvement means that the SHD option can enable much faster, yet reliable characterization of scaled devices at biases where the hot-carrier effects are prominent.

D u c k h y u n C h a n g and Jerry G. F o s s u m

1800

HD Simulation

_.~'.e

.... '.""

:

...'":"-.

~-.

0 < v

t~

o

0

0

SHD Simulation

...-"

....:.,.

-. ,.

...-~" ...-'i...

~

"'",.

8

O

2

< ¢9 O ,.-,. 1

Fig. 5. T w o - d i m e n s i o n a l d i s t r i b u t i o n of the carrier-generation rate in the M O S F E T caused by i m p a c t i o n i z a t i o n as simulated with the H D a n d S H D m o d e l s for Vc,s = VDs = 5.0 V. The predicted s u b s t r a t e currents are 3 4 × 10 ~ A ~ m -~ a n d 3.0 x 10 8 A / ~ m ~ from the S H D a n d H D simulations, respectively.

Simplified energy-balance model for pragmatic multi-dimensional device simulation 40

1801

can be derived as H

~D

dTn'~

2q --~-BB v.'E (9)

dT. = F,(E, T,). dt

(to)

3 {(T, - To) 30

which can be expressed as; 20

lO

o

0.0

012

0,4 016 Voltage Step (&VDs)

0.8

Fig. 6. Average iteration counts for SHD and full HD DC MOSFET simulations done by sweeping VDsfrom 0 to 5.0 V in equal voltage steps A Vos, for different A Vos. As the voltage step increases, the superior computational efficiency of the SHD option becomes prominent.

Now the time-derivative in eqn (10) must be discretized, considering stability and local truncation error in the selection of the numerical integration scheme. First-order backward difference formula (BDF1), second-order backward difference formula (BDF2), trapezoidal rule (TR), and BDF2-TR composite methods have been developed for this purpose[18], and any of these schemes can be used for transient simulation via the SEB/SHD model. For example with BDFI, eqn (10) can be discretized as follows

r ~ - T~-' A~

The SEB/SHD model can be easily extended to a three-dimensional case, for which excessive, usually prohibitive CPU time is required for solving the HD model because of coupling to the third dimension. In the general three-dimensional simulator, the box method, analogous to the two-dimensional one, is commonly used to discretize the carrier transport equations, which are integrated over a small volume enclosing each node[22]. The discretized three-dimensional SEB equation has the same form as eqn (6) with minor changes; length and area in eqn (6) are replaced by area and volume, respectively. This volume can be divided by several triangular prisms. In this case, the current density is defined along each side of each triangular prism, where then the carrier velocity, and thus the energy-relaxation length distributions can be evaluated accordingly from the preceding DD analysis. This capability of extending to the three-dimensional case, potentially with even more enhancement in relative computational efficiency, obviously increases the utility of the SEB/SHD model. So far attention was given to DC or steady-state device simulation. We next discuss extending the SEB/SHD model for transient simulation, and exploiting the benefit afforded with respect to the time discretization. For the general case, this extension is not straightforward; n and p can be important in the time-dependent energy-balance equation, which tends to prevent its simplification. However, the extension can be done effectively for special cases. For example, when T,(dn/ dt)<
(11)

where A~ = t k t k- ~. The use of (11) at each mesh point follows the spatial discretization of eqn (10) as indicated by eqn (6), and the transient SEB/SHD model is then implemented analogously to the DC case as illustrated in Fig. 2. Note that relative to HD, SHD simulation times can be much shorter since larger time steps A~ for transient simulation can be used as implied by the inherent reduced iteration count in Fig. 6. Based on the DC simulations we have done with FLOODS, we project an order-of-magnitude improvement in run time for general two-dimensional transient MOSFET simulation afforded by the SEB/SHD option; and even more for three-dimensional numerical simulation because of the inherent higher degree of nonlinearity in the HD system. -

4. MODEL EXTENSION

- Fn(Ek' Tk)

-

5. SUMMARY

An efficient multi-dimensional SEB/SHD model was developed based on assumed fixed energy-relaxation length distributions, derived from the pre-process DD simulation, in the SEB equations for the SHD simulation. The coupling between electric field and carrier energies is properly taken into account, in contrast to post-processing analytical models and decoupled numerical simulations. However, the coupling between energy and carrier densities, which is typically weak, is implicitly neglected, thereby enhancing the linearity of the system. The new two-dimensional SEB/SHD DC model was implemented in FLOODS as a pragmatic option without substantial source-code modification. All physical models such as impact-ionization and mobility are used without change. Model verification was accomplished through submicron MOSFET simulations. The SHD simulation provides reliable detailed two-dimensional representations for carrier

Duckhyun Chang and Jerry G. Fossum

1802

temperature, carrier velocity, and impact-ionization rate. The accuracy of the model is comparable to the full H D model in F L O O D S . The most noteworthy feature of the new model is in its computational efficiency, which is attributed to the enhanced linearity of the model transport equations. Substantively reduced iteration counts for larger voltage steps, in D C simulation for example, underlies the enhanced efficiency, which in fact can enable pragmatic three-dimensional device simulation. The use of the S E B / S H D model in transient simulation, where quite significant improvements in computational efficiency can result because of large time steps enabled by the enhanced linearity, is not straightforward generally, but follows analogously for special cases. In such cases, and for general DC simulation, the new S E B / S H D option can be easily implemented in any robust two- or three-dimensional numerical H D simulator. work was supported by the Semiconductor Research Corporation. We thank Ming-Yeh Chuang and Mark Law for helpful discussions concerning FLOODS. Acknowledgements--This

REFERENCES

1. Stettler, M. A., Alam, M. A. and Lundstrom, M. S., IEEE Trans. Electron Devices, 1993, 40, 733-740. 2. Scrobohaci, P. G. and Tang, T. W., IEEE Trans. Electron Devices, 1994, 41, 1197 1205. 3. Stratton, R., Phys. Rer., 1982, 126, 2002-2013 4. Blotekjaer, K., IEEE Trans. Electron Devices, 1970, 17, 38-47 5. Tang, T. W., Ramaswamy, S. and Nam, J., 1EEE Trans. Electron Devices, 1993.40, 1469--1477 6. Thoma, R., Emunda, A., Meinerzhagem B., Peifer, H.

J. and Engl, W. L., IEEE Trans. Electron Devices, 1991, 38, 1343-1353 7. Goldsman, N. and Frey, J., IEEE Trans. Electron Devices, 1988, 35, 1524-1529 8. Bordelon, T. J., Wang, X. L., Maziar, C. M. and Tasch, A. F., IEDM Tech. Dig., 1990, pp. 353-356. 9. Chen, D., Kan, E. C., Ravaioli, U., Shu, C. W. and Dutton, R. W., 1EEE Electron Dev. Lett., 1992, 13, 26-28. 10. Bordelon, T. J., Wang, X. L., Maziar, C. M. and Tasch, A. F., Solid-State Electron., 1991, 34(6), 617~28. 11. Crabbe, E. F., Stork, J. M. C., Baccarani, G., Fischetti, M. V. and Laux, S. E., 1EDM Tech. Dig., 1990, pp. 463-466. 12. Slotboom, J. W., Streutker, G., v. Dort, M. J., Woerlee, P. H., Pruijmboom, A. and Gravesteijn, D. J., IEDM Tech. Dig., 1991, pp. 127-130. 13. Souissi, K., Odeh, F., Tang, H. H. K., Gnudi, A. and Lu, P.. IEEE Trans. Electron Devices, 1993, 40, 1501-1507. 14. Agostinelli, V. M., Bordelon, T. J., Wang, X. L., Hasnat, K., Yeap, C., Lemersal, D. B., Tasch, H. F. and Maziar, C. M., IEEE Trans. Electron Devices, 1993, 41, 1784-1795. 15. Cook, R. K. and Frey, J., IEEE Trans. Electron Devices~ 1982, 29, 970-979. 16. Liang. M., Composite thermodynamic system for advanced semiconductor device simulation, Ph.D. dissertation, University of Florida, 1994. 17. Lundstrom, M., Fundamentals of Carrier Transport, Addison--Wesley, New York, 1990 p. 171. 18. Bank, R. E., Coughran, Jr., W. M., Fichtner, W., Rose, D. J. and Smith, R. K., in Process and Device Modeling, ed. W. L. Engl. Elsevier Science, New York, 1986 p. 229. 19, Klaassen, D. B. M., Solid-St. Electron., 1992, 35, 961-967. 20. Ogawa, T., Jpn. J. Appl. Phys., 1965, 4(7), 473-478. 21. Scholl, E. and Quade, W., J. Phys. C: Solid-State Phys., 1987, 20, 861-867. 22. Technology Modeling Associates, DAVINCI User's Manual, TMA, 1992.