ASSESSMENT AND MODELING OF NONQUASISTATIC,
NONLOCAL, AND MULTIDIMENSIONAL EFFECTS IN
ADVANCED BIPOLAR JUNCTION TRANSISTORS
By
JOOHYUN JIN
A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL
OF THE UNIVERSITY OF FLORIDA
IN PARTIAL FULFILLMENT OF THE REQUIREMENTS
FOR THE DEGREE OF DOCTOR OF PHILOSOPHY
UNIVERSITY OF FLORIDA
1992
ACKNOWLEDGEMENTS
I would like to express my sincere gratitude to my
advisor, Jerry G. Fossum, for giving me an opportunity to
work as one of his privileged graduate students on
interesting research topics. Without his devoted guidance,
encouragement, concern, support and patience, this work could
not have reached fruition. My interaction with him has been
a most gratifying learning experience.
I also would like to thank the other members of my
supervisory committee, Professors Dorothea E. Burk, Mark E.
Law, Sheng S. Li, and Timothy J. Anderson, for their
willingness to serve on my committee.
I am also indebted to numerous people I have interacted
with during my stay in Gainesville. First I am grateful to
Mr. D. FitzPatrick for his help in the MMSPICE software
development. Thanks are also extended to many of my
colleagues who helped me through technical discussions or by
cheering me up in difficult times. I cannot mention all of
them, but I should mention Drs. H. Jeong, Y. Kim, J. Choi,
and Messrs. H. J. Cho, S. Lee, H. S. Cho, G. Hong, K. Green,
D. Suh, P. Yeh, M. Liang, D. Apte, S. Krishnan. My deepest
gratitude goes to my parents and sisters Hyesook and Minjung
for their endless love and encouragement throughout the years
of my graduate study. Last but not least, I thank the Lord
for His guidance in my life. I also acknowledge the
financial support of the Semiconductor Research Corporation
and Samsung Semiconductor & Telecommunication Co. Ltd.
iii
TABLE OF CONTENTS
Page
ACKNOWLEDGEMENTS ..........................................
ABSTRACT .. . .............................................
CHAPTERS
1 INTRODUCTION ..........................................
2 MODELING OF MULTIDIMENSIONAL CURRENTS................
2.1 Introduction .................................... ..
2.2 Model Development..................................
2.2.1 Experimental Characterization.................
2.2.2 Analytic Model ................................
2.3 Simulations and Verification ......................
2.4 Summary.......................................... ..
3 NONQUASISTATIC MODELING OF BJT CURRENT CROWDING.....
3.1 Introduction .................................... ..
3.2 Model Development..................................
3.2.1 Switchon Case ................................
3.2.2 Switchoff Case................................
3.3 NQS Model Implementation..........................
3.4 Simulations ......................................
3.5 Summary.......................................... ..
4 ANALYTIC ACCOUNTING FOR CARRIER VELOCITY OVERSHOOT....
.1 Introduction..............................
.2 Model Development........................
4.2.1 Velocity Overshoot....................
4.2.2 Velocity Relaxation...................
4.2.3 Effective Saturated Drift Velocity...
4.2.3.1 Junction SCR......................
4.2.3.2 Currentinduced SCR...............
4.2.3.3 Special case .....................
.3 Comparisons with Energy Transport Model..
.4 Implementation...........................
.5 Simulations ..............................
4 .6 Summary ......................................... 109
5 MMSPICE2 DEVELOPMENT................................. 111
5.1 Introduction ...................................... 111
5.2 New Features.......................................... 112
5.2.1 Multidimensional Currents..................... 112
5.2.2 Current Crowding............................... 113
5.2.3 Velocity Overshoot ............................ 115
5.2.4 Extrinsic Collectorbase Capacitance.......... 115
5.2.5 Substrate Capacitance ......................... 119
5.3 Parameter Evaluation............................... 119
5.4 Model Implementation .............................. 122
5.4.1 Subroutine Modifications....................... 122
5.4.1.1 Subroutine MODCHK.......................... 122
5.4.1.2 Subroutine QBBJT........................... 123
5.4.1.3 Subroutine QBCT ........................... 127
5.4.2 Subroutine Additions .......................... 127
5.4.2.1 Subroutine CROWD......................... 127
5.4.2.2 Subroutine OVERSHOOT ...................... 129
5.5 Demonstration .. .................................. 131
5.6 Summary............................................... 144
6 SUMMARY AND SUGGESTIONS FOR FUTURE WORK................ 145
APPENDICES
A EVALUATION OF JSEO, nEB, JEOP AND nEBP ................... 148
B DISCUSSION ON JQ...................................... 150
C LIMITING JEO(eff) IN THE SWITCHOFF SIMULATION.......... 152
D VALIDITY OF THE DEPLETION APPROXIMATION............... 154
REFERENCES .................. ................................. 156
BIOGRAPHICAL SKETCH........................................ 161
Abstract of Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Doctor of Philosophy
ASSESSMENT AND MODELING OF NONQUASISTATIC,
NONLOCAL, AND MULTIDIMENSIONAL EFFECTS IN
ADVANCED BIPOLAR JUNCTION TRANSISTORS
By
JOOHYUN JIN
August 1992
Chairman: Dr. J. G. Fossum
Major Department: Electrical Engineering
This dissertation is concerned with assessment,
modeling, and simulation of nonquasistatic (NQS), non
local, and multidimensional effects in advanced bipolar
junction transistors. A simple analytic model for the
sidewall injection of the base current, which is shown to be
the most important multidimensional component in scaled
devices, is developed based on the separation of the base
current into internal and peripheral components. Simulation
results for typical test BJTs with various emitter geometries
are compared against corresponding measurements to support
the model. A novel NQS model for transient current crowding
in advanced BJTs is developed for circuit simulation. The
new model, implemented based on a novel use of the previous
timestep solution in the current timestep analysis,
characterizes a timedependent effective bias on the emitter
base junction in a seminumerical analysis, accounting for
base conductivity modulation and the NQS nature of the
crowding. The (dc) debiasing effect, which is important in
analog circuits, is inherently accounted for as well. An
analytic model for electron velocity overshoot resulting from
nonlocal transport in advanced siliconbased BJTs is
developed. The model, which characterizes an effective
saturated drift velocity, larger than the classical value
because of overshoot, is intended for circuit simulation.
The model uses an augmented driftvelocity formalism that
involves a length coefficient derived via Monte Carlo
analysis. The associated velocity relaxation is
characterized phenomenologically to be consistent with
overshoot analysis. The developed chargebased models are
implemented in MMSPICE2, a seminumerical mixedmode
device/circuit simulator, such that users may activate any
combination of the new features by option. The resulting
hierarchical tool, along with the parasitic charge
(capacitance) models included to enhance the usefulness of
the simulator, could indeed enable predictive yet
computationally efficient mixedmode simulations for bipolar
(and BiCMOS) VLSI technology/manufacturing CAD. Utility of
MMSPICE2 is demonstrated by transient simulations of ECL
circuits and devices.
vii
CHAPTER 1
INTRODUCTION
In recent years, advances in process technology have led
to the realization of highperformance bipolar junction
transistors (BJTs). While continual improvement in the
lithographic capability allows the lateral dimensions to be
reduced, scaling down the BJT requires a coordinated change
in both the lateral dimensions and vertical profile to
achieve proper device operation and to improve the intrinsic
device speed. Furthermore, in order to reduce the extrinsic
portion of the bipolar device so that circuit performance can
be more closely tied to the intrinsic device performance,
various selfalignment schemes using polysilicon as base and
emitter contacts have been developed. They all have a
similar structure (see Fig. 1.1), and generally provide much
improved performance over the conventional BJT structure via
a reduction in basecollector junction area and base
resistance.
Despite the impressive progress made in bipolar
technology, computer simulation tools, which are essential to
the optimization of device and circuit designs for the
technology, have not kept pace with it. In integrated
I I Polysilicon
 Nitride
////// Oxide
FisTTT1s +ln*
Fig. 1.1 Cross section of an advanced bipolar junction
transistor fabricated by doublepolysilicon
process.
circuit development and manufacturing today, a technology CAD
(TCAD) system is essential for exploring alternative designs
and evaluating various tradeoffs without timeconsuming and
costly fabrications.
An effective TCAD system requires integrated, physics
based tools for predictive process, device, and (smallscale)
circuit simulation. Computational efficiency is desirable
and indeed essential if the TCAD system is to be used in
manufacturing CAD involving statistical simulation.
Conventional TCAD systems comprise robust, numerical process
and device simulators which drive optimization of empirical
device model parameters for circuit simulation. This
optimization can miss parametric correlations, and hence the
integrated system, although CPUintensive, could yield
nonunique (erroneous) predictions.
Numerical mixedmode device/circuit simulation would
obviate this deficiency, but with a high cost of computation
time. Alternatively, improvement of the TCAD system can
possibly be afforded by incorporation of seminumerical
device models into the circuit simulator which have physical
parameters that relate directly to the device structure. The
resulting tool is an applicationspecific, computationally
efficient mixedmode simulator that can easily be integrated
with the process simulator by a program that evaluates the
model parameters from the doping profile. The MMSPICE
[Jeo90] is such a simulator, which is integrated with SUPREM
[SUP88] by a parameterextraction program, SUMM [Gre90].
The model development for MMSPICE has emphasized the
advanced BJTs. A physical, onedimensional chargebased
model [Jeo89] has been developed and implemented. High
current effects, impact ionization, and nonreciprocal
(trans)capacitances are physically accounted for in the semi
numerical model. This model is sufficient for many
applications, but more work is needed to enhance the
usefulness of MMSPICE.
In most advanced BJTs, the lateral dimension of the
emitter has become the same order of magnitude as the
emitterbase junction depth. Thus, multidimensional current
effects in the peripheral region of the junction are expected
to play a significant role in device performance.
Especially, the variation in commonemitter forward current
gain 3 with geometric shape and size is troublesome to IC
designers [Hwa87]. Hence, some accounting of peripheral
currents is needed for circuit simulation.
Highcurrent effects (e.g., quasisaturation and base
widening, or pushout) are physically accounted for in the
MMSPICE model, but emitter current crowding, caused by
lateral voltage drops in the intrinsic base region, has not
yet been considered. Today's advanced (scaled) BJTs commonly
operate at high current density, and hence transient base
5
current can be much greater than the steadystate current;
this clearly implies the nonquasistatic (NQS) nature of
transient current crowding [Ham88] Therefore, it can be
significant even though dc crowding may be insignificant
[Tan85].
In semiconductor devices where the electric field
increases rapidly over distances comparable to the energy
relaxation mean free path, carrier velocity can overshoot the
value corresponding to the local field because the carrier
(kinetic) energy, which controls the collision time and hence
limits the velocity, lags the field and remains relatively
small [Ruc72]. This nonlocal effect on electron transport
has been recognized as significant in MOSFETs and MESFETs for
years, and now has become important in scaled BJTs [Lee89,
Cra90]. Recent work [Fus92] has indicated that the velocity
overshoot in scaled BJTs can be beneficial, and must be
accounted for in the device and circuit design. However, the
effect has not yet been accounted for in any existing circuit
simulators, and indeed is missing in many device simulators
because of the implied computational intensiveness.
For bipolar integrated circuits, reducing parasitic
capacitances is one of the key issues for speed enhancements.
The extrinsic collectorbase junction capacitance (charge)
has a predominant effect on the circuit performance because
the extrinsic base region is not reduced in proportion as the
intrinsic device is scaled down. The collectorsubstrate
capacitance (charge) is also important.
This dissertation addresses these problems; it is
concerned with the development and implementation of new
models to account for the aforementioned effects in the
advanced BJTs. This work will enable not only truly
predictive, scalable BJT simulations, but also
computationally efficient (seminumerical) mixedmode
device/circuit simulations for bipolar TCAD. The major
contributions made in this work are as follows:
(1) modeling of multidimensional current effects, based on
the separation of the current into internal and
peripheral components;
(2) development of an NQS transient currentcrowding model,
based on a novel use of the previous timestep solution
in the current timestep analysis;
(3) development of an analytic model for electron velocity
overshoot resulting from nonlocal transport in advanced
siliconbased BJTs;
(4) implementation of the new models, including both the
extrinsic collectorbase and collectorsubstrate
capacitances (charges), in MMSPICE to create MMSPICE2.
In Chapter 2, a simple analytic way of accounting for
multidimensional current effects is described. The approach
is based on the separation of the current into areal and
peripheral components. For high VBE, an effective junction
bias (described in Chapter 3) is necessarily defined to
account for the emitter debiasing (a.k.a. crowding) effect.
The model is supported by experimental results of test BJTs
having varied emitter geometries.
In Chapter 3, a new NQS model for transient current
crowding is presented. The model, which characterizes a
timedependent effective bias on the emitterbase junction in
a seminumerical analysis, follows the previous work by
Hauser [Hau64], but physically accounts for base conductivity
modulation and the NQS nature of the crowding. The novel
modeling/implementation is based on the use of the previous
timestep solution in the current timestep analysis, which
in fact could enable general accounting of NQS effects in
seminumerical mixedmode device/circuit simulation. The
tool is supported by numerical simulations of advanced BJTs
using PISCES [PIS84].
In Chapter 4, an analytic model for electron velocity
overshoot in advanced BJTs is presented. The model, which
characterizes an effective saturated drift velocity in the
collector spacecharge regions, is intended for circuit
simulation. The model uses an augmented driftvelocity
formalism that involves a length coefficient derived from
Monte Carlo simulations. The associated relaxation of the
carrier velocity is characterized phenomenologically to be
consistent with the overshoot analysis. Demonstrative
simulation results are presented to assess the significance
of the electron velocity overshoot in advanced bipolar and
BiCMOS technologies, and to support model.
The developed chargebased models are implemented into
MMSPICE2 so that users may activate any combination of the
new features by option. This hierarchical tool is discussed
in Chapter 5. Representative simulations are presented, with
descriptions of the new parameters.
In Chapter 6, the main accomplishments of this
dissertation are summarized, and future research areas are
suggested.
CHAPTER 2
MODELING OF MULTIDIMENSIONAL CURRENTS
2.1 Introduction
For bipolar integrated circuits, reducing parasitic
effects and achieving shallow profiles are two of the key
issues in improving performance. Many selfaligned bipolar
technologies have been developed to achieve low parasitic
capacitance and low base resistivity. They all have a
similar device structure using polysilicon as base and
emitter contacts. In the scaled structure, the distance
between base and emitter contacts is greatly reduced as
determined by the bootshaped sidewall spacer (see Fig. 1.1).
The lateral dimensions of the device have also been scaled
down; for example, the emitter width of today's most advanced
transistors has become the same order of magnitude as the
emitterbase junction depth. Thus, multidimensional effects
in the peripheral region of the junction can play a
significant role in device performance [Hur87].
For digital applications, a most predominant multi
dimensional effect is the lateral injection of significant
base current along the emitter sidewall, which is controlled
by the morphology of the link region [Li88]. One simple way
to reduce this sidewall current component is to increase the
width of the spacer [Dej88, Saw88]. However, many desirable
features of the device depend on the limitation of the spacer
width. For example, as the spacer width increases, the base
resistance and parasitic capacitances increase. Also, the
emittercollector punchthrough current increases due to
insufficient extrinsicintrinsic base overlap in the emitter
periphery [Chu87, Saw88], while an increase in the extrinsic
intrinsic base overlap results in excessive perimeter
tunneling current [Sto83] and hence reduced emitterbase
breakdown voltage. Thus, the control of spacer thickness is
vital to the performance of the device.
The peripheral component of the base current does not
modulate the collector current, and is therefore a parasitic
that degrades the dc current gain 3 in proportion to the
ratio of its magnitude relative to that of the areal
component. Hence, P is degraded more as the perimeterto
area ratio (PE/AE) increases. This implies that the sidewall
effect can be an obstacle for downscaling the emitter size
[Hwa87, Dej88] Therefore, some accounting of peripheral
currents for a given process is needed for a circuit
simulator, e.g., MMSPICE, which actually gives an extra
degree of freedom to the IC designer [Ver87].
In Section 2.2, a simple model based on measurements is
presented to account for the peripheral currents in the
advanced BJT structure. This model, combined with the
currentcrowding analysis described in Chapter 3, will be the
basis for a more predictive and scalable BJT model for
MMSPICE. In Section 2.3, experimental results of test BJTs
having varied emitter geometries are presented to support our
formalism. In fact, interpretation of these results requires
the crowding model of Chapter 3, which was hence developed in
conjunction with the work described in this chapter.
2.2 Model Development
2.2.1 Experimental Characterization
For digital applications, the most important peripheral
current is the sidewall component of the base current.
However, the peripheral component of the collector current is
not significant compared with the areal component, provided
the extrinsic base is welllinked with the intrinsic base
[Li88].
This fact is also supported by our own measurements of
representative (advanced) BJTs provided by Dr. D. Verret of
Texas Instruments. The lateral geometries of the test
devices are described in Table 2.1; LE and WE are the
effective (or actual) length and width of the emitter, and PE
(=2LE+2WE) and AE (=LEWE) are the perimeter and area
TABLE 2.1
LATERAL EMITTER GEOMETRIES OF TEST DEVICES
LE [1Lm] WE [im] PE/AE [1.m1]
9.2 5.2 0.60
9.2 4.2 0.69
9.2 3.2 0.84
9.2 1.7 1.39
9.2 1.2 1.88
9.2 0.7 3.08
9.2 0.45 4.66
respectively. The spacer width of these devices is estimated
to be 0.4j.m. Fig. 2.1 shows the base (JB) and collector (Jc)
current densities versus PE/AE for the devices with LE fixed
at LE=9.2RLm when VBE=0.4 or 0.7V. Since JC is almost constant
regardless of PE/AE as well as VBE and VBC, we infer that the
peripheral collector current can be neglected at least for
relatively low VBE. On the contrary, JB clearly increases
with PE/AE, obviously implying a significant lateralinjection
component. We note that this parasitic effect becomes more
significant as VBE is reduced, which we believe reveals that
the peripheral base current is due to the recombination of
excess carriers in the peripheral junction spacecharge
region (SCR) near or at the oxidesilicon interface.
The lateral injection can be understood better if the
peripheral component of base current is quantified.
Empirically, the total base current IB can be separated into
areal and peripheral parts as follows [Rei84]:
IB = IBA + IBP
_VBE VBEv
SCAA e xp 1 + CPPE exp 1 (2.1)
where CA, nA, Cp, and np are (processdependent) empirical
constants, which can easily be evaluated using the basic
experimental method discussed in Appendix A. In (2.1), the
E
0
"l
13
Oa
1010
1011
10
U j I I I I I I I I I I I I I III I
0 1 2 3 4 5
PE/AE [1/um]
(b)
E
m
o
rn
106
107
1087
18
PE/AE [1/um]
Fig. 2.1 Base and collector current densities versus PE/AE
for devices with LE=9.2tm: (a) VBE=O .4V; (b)
VBE=0.7V.
6 6 6 0
Sa
A VBC= 3.0V
O VBC= O.OV
6 6 6 6
S VBC = 3.0V
0 VBC = O.OV
0
voltage drop across the extrinsic base resistance is
neglected for lowcurrent conditions.
Based on this formalism, it is possible to calculate the
contribution of the peripheral current to the total base
current. Doing this for the devices previously characterized
yields in Fig. 2.2 IBP/IB versus PE/AE for VBC=O.OV. As
discussed before, the peripheral base component increases
with PE/AE. For example, when VBE=0.7V and PE/AE=0.60/Llm
(actually, this is equivalent to the device with WE=5.21m),
IBP is only 16% of the total base current, but it increases to
50% when PE/AE=3.1/tm (i.e., WE=0.7 m). For reduced VBE, the
effect of lateral injection becomes more significant in
accord with our previously stated recognition; when VBE=0.4V,
the mentioned ratios are changed to 47% and 82% respectively.
Our other simulations and measurements show that the
peripheral collector current evaluated via this methodology
is about 10% of the total collector current on the average.
2.2.2 Analytic Model
With this insight, we can extend the MMSPICE BJT model
to account for the peripheral base region, at least to first
order. The extended model is restricted to include only the
lateral injection of the base current, which has been shown
to be the most important multidimensional effect in modeling
I I 1 1 1I I 1 1I I 1I r 11 11
1
0.8
0.6
0.4
0.2
0
1 2 3 4 5
PE/AE [1/um]
Fig. 2.2 Simulated IBP/IB versus PE/AE for the devices used
in Fig. 2.1.
VBE =0.7V
BE
advanced BJTs. Based on the insight derived from the
measurements, we add only a peripheral component of base
current to the existing BJT routine in MMSPICE. This
additional component is proportional to the emitter perimeter
PE, and represents peripheral SCR recombination near the
surface. The peripheral base current IBp can be expressed as
IBP = JEOPPE exp VBE 1 (2.2)
L \nEBpVT
where JEOP and nEBP represent the peripheral saturation
current density (per unit length) and the peripheral emission
coefficient respectively. The sidewall injection effect
could also be dependent on the emitter junction depth, but we
assume that this dependence is implicitly included in the
above formalism.
In a dc case, the predominant components of the areal
base current are typically backinjection current from the
base to the emitter and the recombination current at the
(emitterbase) junction SCR. (Recombination in the quasi
neutral base and the epi collector is neglected here since it
is typically insignificant in advanced BJTs.) Hence, the
total base current IB can be expressed as
IB = IBA + IBP
= JEoAE [expB 1 + JSEOAE [exp( VBE
T nEBV
+ JEOPPE [exp V1E (2.3)
L \nEBpVT J
where JEO is the (areal) emitter saturation current density,
and JSEO and nEB are (areal) SCR saturation current density
and SCR emission coefficient respectively.
Although (2.3) is sufficient for many operating ranges,
it is necessary to examine whether it is valid for high
current operation where additional effects are significant.
In this case the actual (peripheral) junction bias V'BE cannot
be approximated as the terminal voltage VBE; V'BE is
considerably less than VBE since the voltage drops across the
extrinsic base and emitter resistances are no longer
negligible. Furthermore, the areal component is degraded by
the lateral voltage drops in the intrinsic base region. In
fact, interpretation of data necessitated the current
crowding modeling described in Chapter 3. Hence we modify
(2.3):
IB = JEOAE expBE(e) 1 + JSEOAE [exp(V ) 
VT nEBV E
+ JEOPPE expVIB E 1 (2.4)
L \ngBpVTl I
where VBE(eff) is defined (in Chapter 3) as the effective bias
on the emitterbase junction to account for the debiasing
(a.k.a. current crowding) in terms of the actual (peripheral)
bias V'BE. Note that in (2.4), the peripheral current term is
not threatened by the current crowding because the peripheral
junction voltage is always fixed at V'BE. Although the
debiasing effect was classically characterized by Hauser
[Hau64], his treatment is inadequate for advanced BJTs
because it neglects conductivity modulation of the base. On
the contrary, the concept of the effective bias can account
for the highcurrent effects via the chargebased BJT model
[Jeo89]. When the debiasing effect is significant, the
effective bias is of course less than the actual junction
bias V'BE. (In this case, V'BE is also significantly less
than VBE.) Otherwise, VBE(eff) would be almost the same as
V'BE. This effective bias is derived from the quasithree
dimensional crowding analysis, which involves a coupling of
the vertical and lateral carriertransport analyses in the
base region. Details are described in Chapter 3.
Fig. 2.3 illustrates (V'BEVBE(eff))/VT versus WE
predicted by the debiasing analysis for typical advanced
devices with LE=9.211m. When VBE=0.7V, the debiasing effect
is, as expected, negligible resulting in VBE(eff)=V'BE=VBE
regardless of WE and VBc. However it becomes noticeable for
higher VBE and especially for greater WE, due to the increased
0.8
0.6
 
S0.4 
w
m
0
0.2
I t H I ~ ~ I I I I i I I I I I 1 1
VBC=3.0V
VBC=O.OV
VBC=O.O or 3.0V
.. .. ..
SI . . . . I I
0 1 2 3 4 5 6
WE [urn]
Fig. 2.3 Simulated (V'BEVBE(eff))/VT versus WE for typical
advanced BJTs with LE=9.21m.
I . I I I I I I i I , I I I I I I I I I I I ,
voltage drops in the intrinsic base region. The debiasing
effect also becomes more important with increasing reverse
bias on the basecollector junction because the base
resistivity increases correspondingly. For contemporary
scaled BJTs however, it is not significant [Tan85]; for
WE=2jm at VBE=0.9V and VBC=3.0V, the voltage difference
between the actual and effective bias is about 20% of the
thermal voltage.
2.3 Simulations and Verification
The test devices, representative of the advanced bipolar
technology, were used to verify the model. The devices, from
Texas Instruments, were fabricated using a doublepolysilicon
process in conjunction with a sidewall spacer technique,
which enables a selfaligned submicrometer emitter structure.
In order to identify significant multidimensional effects,
transistors with different PE/AE (see Table 2.1) were
measured.
Simulations were done with MMSPICE2, which includes the
peripheral base current [eq. (2.2)] and the currentcrowding
model as described in Chapter 3. At first, the model
parameters associated with the lateral injection were
extracted as described in Appendix A. Then, with no
additional parameter extraction, all BJTs were simulated with
reasonably good accuracy simply by scaling AE.
Simulated Ic/WE and IB/WE compare quite well with the
corresponding measurements in Fig. 2.4(a) when VBE=0.4V and
VBC=O.OV. Note that the lateral injection effect on the base
current becomes significant as WE is scaled down; IB/WE
increases because the ratio of the peripheral to the areal
component increases. However, the contribution of the
peripheral collector current is negligible for each device.
Note that if IBP had not been accounted for, IB/WE would have
been predicted to be a constant, since the voltage drops
across the extrinsic resistances are negligible for each
device at this bias point. For the corresponding 0 shown in
Fig. 2.4(b), the simulations are excellent. As expected, P
is reduced with decreasing WE. Although 3degradation is an
obstacle for downscaling WE, we expect that our firstorder
accounting of the lateral injection could give an extra
degree of freedom to the circuit designer.
The peripheral collector current is still negligible
when VBE is increased to 0.7V, as shown in Fig. 2.5(a).
Still, the sidewall injection of the base current, although
not as significant as in the lowcurrent region, is important
especially for devices with small WE. The simulations are
good, although there is a small discrepancy between the
measured and predicted values of IB/WE for submicron devices.
Indeed this discrepancy seems to be inevitable because the
10100
LUi
S1011
10 11
140
120
100
an
60
40
20
0 1 2 3 4 5 6
WE [um]
(b)
L''~' , I , I , I , I I ,
A VBC= 3.0V
O VBC= 0.0V
 Simulation
1 2 3 4 5
WE [um]
Fig. 2.4 Measured and simulated IC/WE, IB/WE in (a) and P in
(b) for the test BJTs with LE=9.2.m for VBE=0.4V.
O Measurement (VBC=O.OV)
Simulation
00n o n 0
0
(a)
Ir1 1, 1 1, 1
105
106
107
108
0
O
0
~ VBC= 3.0V
0 VBC= O.OV
 Simulation
0 1 2 3 4 5 6
WE [um]
Fig. 2.5 Measured and simulated IC/WE, IB/WE in (a) and P in
(b) for the test BJTs with LE=9.2pm for VBE=0.7V.
2 3 4
WE [um]
Oo o o n O
0 Measurement (VBC=O.OV)
Simulation
0
I
r r a I I I I i s r r r I
f i l l
IBprelated model parameters were evaluated from the devices
operating in the low current region; according to (2.3), the
PEdependent term would become negligible with increasing VBE.
However, our model seems adequate, as implied by the
corresponding 3 results in Fig. 2.5(b).
For VBE=0.9V in Fig. 2.6(a), the simulations are also
reasonably good. We note that IC/WE and IB/WE decrease with
increasing WE, not because the lateral injection becomes less
significant as in Figs. 2.4 and 2.5, but because both the
debiasing of the internal junction and highcurrentinduced
voltage drops across the extrinsic resistances, including
base resistance, increase with WE. From the figure however,
we can infer that the voltage drops, which become greater for
large devices due to the increased terminal currents, are
most dominant. The effect of current crowding on P is well
illustrated in Fig. 2.6(b); of course, the better simulations
obtain with debiasing accounted for. However the debiasing
seems to be insignificant for contemporary scaled devices, as
discussed before. Our other simulations show that for
devices with WE>LE, the debiasing effect is almost the same
for each device, since the predominant base current flow
under the rectangular emitter is laterally along the shorter
emitter dimension (LE in this case).
103
0 1 2 3 4 5 6
WE [um]
(b)
. . I . . eI I I I * si, ,
140
120
100
80
60
40
20
2 3 4
5 6
WE [um]
Fig. 2.6 Measured and simulated IC/WE, IB/WE in (a) and 3 in
(b) for the test BJTs with LE=9.2lm for VBE=0.9V.
0
 
O Measurement (VBC=0.(
  w/o Crowding
S w/ Crowding
0
VBC = 3.0V
O VBC = O.OV
 w/ Crowding
   w/o Crowding
6p
2.4 Summary
A simple analytic model for the lateral injection of
base current, which is shown to be the most predominant
multidimensional current effect in advanced BJTs, has been
developed by separating the base current into internal and
peripheral components. The model is intended for (digital)
circuit simulation and has been implemented in MMSPICE2.
For high VBE, the effective bias (see Chapter 3) on the
emitterbase junction is defined to account for the debiasing
effect. The tool is well supported by experimental results
of test BJTs having varied emitter geometries, despite the
fact that the simulation for each device was done by scaling
only AE for a given parameter set. Therefore, this lateral
injection model, combined with the currentcrowding analysis,
can be the basis for more predictive and scalable BJT
simulation for TCAD.
For analog circuit simulations, more precision is
usually required. In this case, it is possible to analyze
more physically the multidimensional effects by cascading a
second (peripheral) BJT to the intrinsic one, each
represented by the onedimensional BJT model in MMSPICE; the
composite transistor is also useful to account for the
parasitics associated with the extrinsic base region as well
as the lateral injection effect, for example in RF IC design
applications [Jaf92].
From the measurement and simulation results for
contemporary BJTs, the following conclusions were reached:
(1) The lateral injection of the base current becomes more
significant with decreasing VBE, which reveals that the
nature of this perimeter effect is recombination at the
peripheral junction SCR near the oxidesilicon interface.
(2) The peripheral component of collector current is
typically negligible.
(3) In highcurrent regions, the voltage drops across the
extrinsic resistances are most predominant, and the dc
debiasing effect seems to be negligible for contemporary
BJTs.
CHAPTER 3
NONQUASISTATIC MODELING OF BJT CURRENT CROWDING
3.1 Introduction
In contemporary digital circuits containing advanced
(scaled) BJTs, high transient base current can be much
greater than the steadystate current; this clearly implies
the nonquasistatic (NQS) nature of transient current
crowding. (We generally define an NQS effect in the time [or
acfrequency] domain as one that cannot be inferred nor
characterized from steadystate [dc] conditions.) Hence it
can be significant even though dc crowding may be
insignificant [Tan85]. The classical treatment of emitter
current crowding by Hauser [Hau64], although useful, is
inadequate for advanced BJTs because it neglects conductivity
modulation of the base, which can occur because of high
injection and/or base widening, and because it assumes
steadystate or quasistatic conditions. In fact, transient
current crowding is NQS, as well as being dependent on the
base conductivity modulation [Ham88].
There has been some modeling done addressing the NQS
nature of current crowding, but generally involving
distributed lumpedmodel representations of the base region.
Indeed NQS effects can be physically accounted for by
cascading a sufficient number of elemental quasistatic
models, but computational efficiency must be sacrificed. Rey
[Rey69] used a more novel approach to model ac crowding and
derived a frequencydependent base impedance for an
equivalentcircuit model.
In this chapter we extend the onedimensional BJT model
in MMSPICE1 to account for threedimensional transient
current crowding in advanced, selfaligned devices which have
peripheral base contacts. The formalism includes a novel
methodology for seminumerically modeling general NQS effects
in transient device/circuit simulation. The new model
characterizes a timedependent effective bias on the emitter
base junction for each NewtonRaphson iteration of the
circuit nodal analysis at each timestep. The seminumerical
analysis follows Hauser, but physically accounts for base
conductivity modulation and the NQS nature of the crowding.
The latter extension is effected by the novel
modeling/implementation that involves the use of the previous
timestep solution in the current timestep analysis. The
model naturally accounts for dc crowding as well, which is
important in analog circuits, and which was needed in Chapter
2 to interpret the multidimensional current measurements in
the BJT. It does not require a lumped intrinsic base
resistance [Jo90], which is commonly used in BJT circuit
models.
The NQS model, implemented in MMSPICE2, enables a semi
numerical mixedmode device/circuit simulation capability for
applicationspecific TCAD. The tool is supported by
numerical simulations of advanced BJT structures using PISCES
[PIS84]. It is used to clarify the nature of the added (NQS)
delay due to current crowding in switchon and switchoff
transients in representative BJT inverting circuits, and it
reveals the significance of transient crowding even in
submicron devices.
3.2 Model Development
The intrinsic base of the advanced (selfaligned) BJT is
surrounded by a highconductivity extrinsic base. Hence the
predominant base current flow under a rectangular emitter is
along the shorter emitter dimension (WE); this is assumed in
our (quasithreedimensional) crowding analysis. Consider a
section of the base of an npn BJT as shown in Fig. 3.1, where
WE is shorter than the emitter length LE. For transient
conditions at a point in time, let iB(y) be the lateral base
current which causes the crowding in the emitterbase
junction. Then, the emitterbase junction voltage v(y) can
be expressed as
EMITTER
I Y
y=O
y=WE/2
Fig. 3.1
Cross section of the advanced (symmetrical)
bipolar junction transistor. Wb(eff) is the
widened (due to possible quasisaturation) base
width.
BASE
v(y) = v(O) dv
= VBE iB(y)dRBi
= VBE iB(y)pdy (3.1)
where VBE is the peripheral junction voltage and p is the
specific base resistivity,
dRBi 1 WE
 (3.2)
dy 2pqLpLEWb(eff) 2jLp(QBB + QQNR)
In (3.2), p represents an average hole density at y, which we
assume can be represented in terms of the total hole charge
(QBB+QQNR) in the quasineutral base (possibly widened to
Wb(eff) due to quasisaturation); QBB, the hole charge in the
metallurgical base region, and QQNR, the hole charge in the
widened base region, both integrated over the emitter area AE
as well as over the base width, are characterized in the one
dimensional model [Jeo89] This assumption in (3.2) is
consistent with a quasitwodimensional analysis (to be
described) which links the onedimensional ambipolar
transport to the lateral hole flow. Implicit in the
assumption is a neglect of lateral hole diffusion, which
indeed is typically small compared to the lateral drift
current when crowding is significant. The model deficiency
resulting from this neglect will be shown to be
inconsequential later. The hole mobility at y is also
approximated by an average value Lp which is reasonably
estimated from common sources. Note that the factor of 2 in
the denominator of (3.2) accounts for the symmetry of the
transistor obvious in Fig. 3.1.
For transient excitation, the main components of the
intrinsic base current iB are typically hole current back
injected from the base to the emitter (IBE) and majorityhole
charging/discharging current (dQBE/dt) Note that QBE
includes components of (QBB+QQNR) communicating with the
emitter [Jeo89]. It comprises space charge (e.g., junction
depletion charge) as well as quasineutralregion charge in
the intrinsic device structure. Generally, IBE(y) can be
expressed as
IBE(y) = IBE(O) 2JEoLEexp ) l dy (3.3)
where JEO is the (constant) emitter saturation current
density. We assume that the ydependence of dQBE(y)/dt, at a
particular point in time, may be similarly expressed as
m+l m+l1
dQBE (y) dQBE (0) j v(y)
BE = B 2J L exp v ( 1dy
dt dt f Q VT 1 (3.4)
where JQ(t) is a transient (timedependent) counterpart to
JEO. Implicit in (3.4) is an idea that JQ can be estimated
from the previous timestep (t=tm) solution for dQBE/dt for
use in the current timestep (t=tm+l) analysis as follows:
dQmE (0)
Ja+l = dt
SLEWf BE(eff) (3.5)
VT
where VBE(eff) is an NQS effective bias on the emitterbase
junction defined (see (3.7)) to account for the current
crowding (see the discussion in the Appendix B). So our
model, when implemented based on the previous timestep
solution, accounts for transient crowding nonquasi
statically. The approximation in (3.5) is viable even for
fast transients because of the automatically controlled time
step reduction in the simulator, which is needed to ensure
acceptable truncation error and convergence of the timepoint
solution.
With (3.3) and (3.4), the intrinsic base current iB(y)
is written as
dQBE(y)
iB(y) = IBE(Y) + dQBE
dt
d (0y
= IBE(0) + dQBE(0 2(JEo + JLE expP 1 dy
dt f VT
= B(O) 2JEO(eff)L exPv 1dy (3.6)
where the timedependent JEO(eff) is defined as the sum of JEO
and JQ(t) To facilitate an analytic accounting for the
crowding (reflected by the integral in (3.6)), we define
VBE(eff) based on the total intrinsic base current:
VBE (eff) 1
iB(0) = LEWEJEO(eff) exp V 1 (3.7)
Note that (3.7) is consistent with (3.5).
Now, following Hauser's classical analysis [Hau64], we
differentiate (3.6) combined with (3.1) to get
S= 2JEO(eff)L exp(VBE iB(y)pdy) 1 (3.8)
Dy VT
This integraldifferential equation for iB (y) may be
transformed into a closedform secondorder differential
equation by differentiating it. This differentiation, with
v(y)
ex >>[ 1 (3.9)
VT
for all values of y, which is generally valid for problems of
interest, yields
2 iB P iB
+ B (3.10)
ay2 VT ay
For transient crowding, (3.10) has two different types of
solution depending on the sign of BiB/ay. We consider the two
cases separately.
3.2.1 Switchon Case
When the BJT is switchedon, iB>O tends to cause
peripheralemitter current crowding, as in dc crowding
[Tan85]. In this case, aiB/By is negative, and the solution
of (3.10) is
iB(y) = A tan _B B (3.11)
L2V \ B j
where A and B are arbitrary constants of integration. The
constants can be evaluated from the boundary conditions of
the problem. For the structure shown in Fig. 3.1, we have
due to the symmetry
iWE ,
X21
which gives B=WE/2. Then from (3.11),
iB(y) = A taz( 2y
W \
(3.12)
(3.13)
where z=ApWE/(4VT) Hence, the total base current is
iB(0) = A tan(z)
(3.14)
which is equated to (3.7) to characterize vBE(eff)
Using (3.13) in (3.1) and doing the integration yields
v(y) = VBE 2VT ncsz(l 2y/WE)
cosz
(3.15)
Note for this case that
v(W) = VBE + 2VT In(cosz) < VBE = v(0)
Now using (3.15) in the integration in (3.6), with the
boundary condition (3.12), yields another expression for the
total base current:
BE sinz COSZ
iB(O) = LEWEJEO(eff)exp  sin (3.16)
With (3.7),(3.14), and (3.16) we now have a set of three
nonlinear equations in three unknowns (vBE(eff) A, and
iB(O)), which can be numerically solved by the iterative
NewtonRaphson method. An interesting relationship is an
expression relating VBE(eff) to VBE. This is obtained by
equating (3.7) to (3.16):
xVBE(ef f)) V=BE sinz COSZ
exp = exp (3.17)
VT VZ z
for exp [BE(eff)/VT]>>i. Note that vBE(eff) is always less
than VBE for the switchon case since (sinz cosz/z) is less
than unity.
The accounting for dc crowding in the model is inherent
in the switchon analysis described above. For the dc case,
JQ=O and JEO(eff)JEO in (3.16)
3.2.2 Switchoff Case
For the switchoff case, iB<0 tends to cause central
emitter current crowding. The analysis is very similar to
that for switchon, except that now aiB/ay is positive.
Actually this condition does not obtain instantaneously when
the BJT is abruptly turned off from an onstate. A very fast
transient occurs during which holes diffuse out of the
intrinsic base periphery to support the centralemitter
crowding that ultimately controls the predominant switchoff
transient. Our model presented below is invalid during this
fast transient since it neglects lateral diffusion flow.
However this brief invalidity is typically inconsequential
with regard to simulating the predominant transient. Note
that the fast (diffusion) transient is governed by a lateral
quasineutral base transit time for minority electrons; it is
proportional to (WE/2)2/Dn where Dn is an average diffusion
constant for electrons.
With the same boundary condition (3.12), the solution of
(3.10) with iiB/ay>0 is
iB(y) = A tanhz1 2y (3.18)
So, the total base current is now
iB(O) = A tanh(z) (3.19)
Once again we define the NQS effective bias vBE(eff) by (3.7),
in which JEO(eff) is now negative because predominant
discharging current flows in this case. Following the steps
in the switchon analysis, we get
v(y) = VBE + 2VT in[cosh{z( S 2y/z (3.20)
cosh{z(1 2y/WE)
Note here that
v = VBE + 2VT In(coshz) > VBE = v(0)
The total base current can now be derived, analogously to
(3.16), as
VBE COshz sinhz
iB(O) = LEWEJEO(eff)exp(T coshz sinhz
T (3.21)
Once again we have a system of three nonlinear
equations, (3.7), (3.19), and (3.21), that define VBE(eff)
seminumerically via iterative solution. Another interesting
relationship between VBE and vBE(eff) is obtained from (3.7)
and (3.21):
VxBE(eff) exVBEi coshz sinhz
exp = exp (3.22)
VT VT z
for exp[vBE(eff)/VT>>1 Note that vBE(eff) is always greater
than VBE in the switchoff case since (coshz sinhz/z) is
greater than unity.
We note that the switchoff analysis described above has
no solution for extremely large negative JEO(eff), which tends
to obtain when the discharging current dQBE/dt (viz., JQ in
(3.5)) becomes too large compared with the dc current IBE*
This condition is nonphysical, and reflects the deficiency
of our model during the initial fast (diffusion) transient
discussed previously.
The nosolution problem can be avoided by limiting
JEO(eff). Such limitation results in a solution, albeit
invalid, that most importantly carries the simulation through
the fast transient to the most significant lateraldrift
controlled switchoff transient. So, for each iteration at
each timestep, we calculate a hypothetical maximum absolute
value of JEO(eff) for which the system of equations is
solvable, and then compare it with the actual JEO(eff); the
smaller value is used for the analysis. Details are given in
the Appendix C. This hypothetical limit for JEO(eff) is, as
expected, used only at the very beginning of the switchoff
transient, where the model is nonphysical anyway, and indeed
is insignificant with regard to the predominant transient.
3.3 NOS Model Implementation
Our novel NQS modeling/implementation in MMSPICE2 of
the BJT current crowding involves a coupling of the vertical
and lateral carriertransport analyses in the base region.
For the npn device, the analysis of the twodimensional hole
flow seminumerically defines vBE(eff) for each NewtonRaphson
iteration of the circuit nodal analysis at each time step.
The implemented transientcrowding model algorithm is
flowcharted in Fig. 3.2. The calculation of JQ from the
previous timestep solution for use in the current timestep
is done only in the first iteration at each time step, and
the value is used for all subsequent iterations. With the
terminal biases VBE and VBC passed in from the nodal analysis,
the onedimensional model routine in MMSPICE solves the
ambipolar transport, accounting for constant extrinsic
terminal resistances, and characterizes the base charge in
both the metallurgical (QBB) and widened (QQNR) base regions.
These charges define the specific base resistivity (p) for the
current timestep analysis, which is needed in the solution
of the hole transport to derive a new VBE(eff). As discussed
in Section 3.2, this derivation requires a NewtonRaphson
iterative solution because of nonlinearities due to the
conductivity modulation.
Note in Fig. 3.2 that vBE(eff) is not iteratively coupled
to the onedimensional model solution; that is, p is not
updated to correspond with vBE(eff) VBE. Although this one
pass derivation of vBE(eff) using p(VBE) might seem incomplete,
it is proper. A complete iterative solution, which would
Fig. 3.2
Flowchart of the MMSPICEimplemented transient
current crowding analysis, for every iteration at
each time step.
require an outer Newtonlike loop in the algorithm, would be
nonphysical. The reason is that in the switchon case where
VBE(eff) is less than vBE, the smaller VBE(eff) in the one
dimensional model would not adequately account for possible
highcurrent effects at the periphery, and that in the
switchoff case where VBE(eff) is greater than VBE, the larger
VBE(eff) in the onedimensional model would tend to diminish
the central crowding effects by implying a smaller p.
With VBE(eff), the onedimensional MMSPICE model routine
is called again to obtain the nominal biaspoint solution.
Since the model is seminumerical, analytic derivatives of
the currents and charges cannot be given explicitly. Thus,
numerical (divideddifference) approximations are used to
evaluate (trans)conductances and (trans)capacitances for
use in the subsequent nodal analysis. In order to do that,
the model routine is called twice more with perturbed values
of vBE(eff) and vBC as indicated in Fig. 3.2. The admittance
matrix is then loaded, and ordinary circuit nodal analysis
follows.
3.4 Simulations
Examples of transient simulations using MMSPICE2 are
presented in this section. One circuit chosen for simulation
is a singletransistor inverter shown in Fig. 3.3, with no
Vcc=2 V
RCC=200 Q
0 OUT
RBB=100
Fig. 3.3
A single transistor inverter circuit. The base
terminal is driven with a voltage pulse that is
delayed by 200ps and then ramped up (down) from
0.4V (0.9V) to 0.9V (0.4V) at a rate of 0.1V/ps.
load on the output. The assumed BJT model parameters
characterize a typical advanced device structure with
WE=1.2tm. The peak base doping density is 1.5x1018cm3 and
the metallurgical base width is 0.15p.m. For the switchon
transient, the NQS nature of the transient current crowding
is well illustrated in Fig. 3.4 where the simulated time
dependent JQ, defined in (3.5), is compared with JEO. Note
that JQ is several orders of magnitude greater than JEO at the
moment the device is switchedon. It decreases monotonically
with time and finally becomes less than JEO only when the
device nears steady state. In the switchoff case, JQ is
negative, and its magnitude is not so large as for the
switchon case. This is due to the exp[vBE(eff)/VT] term in
the denominator of (3.5), which is large when the device is
switched off.
For the complete switchon/switchoff cycle, Fig. 3.5
contrasts the simulated vBE(eff) with VBE in time, accounting
for constant extrinsic/external base resistance, which is
reflected by the discrepancies between VBE and the input
voltage vin. The moment the device is switchedon, vBE(eff)
becomes, as mentioned earlier, less than VBE due to the high
transient base currentinduced crowding, but then increases
steadily with time to a value that corresponds to dc
crowding, which is relatively insignificant. For the switch
off transient, vBE(eff) is greater than vBE, but the difference
48
I .I , I I , I I I
0 11010
21010
3 1010 41010
Time [sec]
51I 610 I 710
51010 61010 71010
Fig. 3.4
Simulated JQ versus time in the switchon case.
JEO is the emitter saturation current density.
102
10
10 
2
10
i
106
108
0
0
0
o
o 
o
o
J 00o
  E O 0 0 0 00 0 0
O000
4 .
. . .
. . l ' '
1.2
1
LU 0.8
oa
" 0.6
LU
0.4
0.2
0 21010 41010 61010 81010
1 109
Time [sec]
Fig. 3.5
Simulated VBE(eff) versus time for the complete
switchon/switchoff cycle. The input pulse and
the actual (peripheral) baseemitter junction
voltage are shown for comparison.
is not so noticeable as for the switchon case. These
results suggest that the centralemitter current crowding
during a switchoff transient is much less significant than
the peripheralemitter crowding during a switchon transient.
This can be attributed to the level of base conductivity
modulation (reflected by p) at the initial stages of the
respective transients.
Fig. 3.6 shows the output voltage characteristics of the
inverter simulated with (MMSPICE2) and without (MMSPICE1)
the current crowding accounted for. In accord with
conclusions drawn from Fig. 3.5, the result of the switchon
transient crowding is a substantively slower response, while
the added delay is insignificant for the switchoff
transient. Other simulations show that accounting for only
quasistatic crowding (due to JEO in (3.6)) yields an output
voltage characteristic which is virtually identical to that
predicted by the simulation in Fig. 3.6 for which crowding
was completely neglected.
Predicted switchon delays of the single transistor
inverter versus WE, with the emitter area fixed
(AE=LExWE=9.2x2.0m2), and with the emitter area scaled with
WE, are plotted in Fig. 3.7. The emitter width WE was varied
using the values 0.1, 0.4, 1.2 and 2.Opm. The delay was
defined as the time for the output current to reach 50% of
its final (high) value. The effect of the crowding is made
21010 41010 61010 81010
Time [sec]
Fig. 3.6
Output voltage characteristics of the single
transistor inverter simulated with (MMSPICE2) and
without (MMSPICE1) the transient current crowding
accounted for.
1 109
I I i I i I
' ' I '.. I .1 I I ' ' I
0 0.5 1 1.5 2 2.E
WE [um]
Fig. 3.7
Predicted switchon delays of the single
transistor inverter versus WE, with the emitter
area fixed (AE=9.2x2.0lm2), and with the emitter
area scaled with WE.
MMSPICE2, AE fixed
  MMSPICE1,AE fixed
[ MMSPICE2, AE scaled
 f MMSPICE1, AE scaled
200
160
0 120

a 80 
40
0
.... I
 )  O
apparent by including in the figure delays predicted by one
dimensional (MMSPICE1) simulations. For the switchon
transient, the results, consistent with previous work
[Tan85], show that peripheralemitter crowding causes an
added delay, one that tends to become insignificant only when
WE is reduced to deepsubmicron values [Ham88]. Note in Fig.
3.7 that when the emitter area is scaled with WE, the delay
is more sensitive to WE. The reason of course is that, in
addition to the crowding effect, the amount of charge that
must be stored in the BJT varies with WE. Other simulations
show that the relative importance of the crowding varies
inversely with the extrinsic (plus external) base resistance.
Results of switchoff simulations with varying WE show
that the added delay due to centralemitter crowding is
negligible, at least for WE<2p~m. Indeed the simulations
predict that the reduced delay of a scaled (WE and AE) device
is due predominantly to the reduced charge storage in the
BJT.
The effect of the emitter length LE on the current
crowding is reflected in Fig. 3.8, which shows normalized
predicted switchon delays versus WE for devices with AE fixed
at 9.2x2.Om2 or 3.2x2.Om2. Note that for a fixed WE, the
crowding effect on the delay diminishes with increasing LE.
This is due to the decreasing specific resistivity p in
(3.2).
0 0.5 1 1.5 2 2.5
WE [um]
Fig. 3.8
Predicted normalized switchon delays versus WE,
with fixed AE, for devices with different LE.
The influence of the nominal base resistivity, viz., the
Gummel number, on the added switchon delay due to crowding
is revealed in Fig. 3.9 where predicted normalized delays are
plotted versus WE (with fixed AE=9.2x2 0m2) for three
different metallurgical base widths WBM. The peak base doping
density was fixed at 1.5xl018cm3. The plots show how the
transient crowding becomes more significant as WBM is scaled
down, independent of the increasing current gain of the BJT
since there is no load on the inverter (Fig. 3.3).
In order to verify our model, twodimensional numerical
simulations of the nominal BJT inverter were performed using
PISCES [PIS84], the results of which for varying WE are shown
in Fig. 3.10. In these switchon and switchoff simulations,
the actual emitter length was fixed at lLm because the output
currents of PISCES are always normalized by the length
perpendicular to the simulated structure. Also, the values
of WE used for the plots are the effective emitter widths,
which are about 0.2pLm wider than the polyemitter windows
because of lateral diffusion. The contact resistances at the
collector and base terminals were specified to include the
external resistances in the inverter circuit. Included in
Fig. 3.10 are corresponding MMSPICE device/circuit
simulations, with LE=1lm. In the switchon case, the
transient current crowding is significant and is faithfully
predicted by MMSPICE2, as contrasted by the inaccurate
,, ,I I ~ I I I _
S1 11 1 1 I I 1 I 1
0 0.5 1 1.5 2 2.5
WE [um]
Fig. 3.9
Predicted normalized switchon delays versus WE,
with fixed AE, for devices with different WBM.
1 i I I I I I I
1.2
1 
0.8
0.6
0.4
0.2
 WBM=0.20um
  WBM=0.15um
WBM=0.10um
" '
' '
140
120 
100
80
60 
20 
0
0 I
0.5
WE [um]
Fig. 3.10
PISCES simulations of the switchon and switchoff
delays versus (effective) WE of the single
transistor inverter, with corresponding MMSPICE
simulations. LE=1~m for all simulations.
0 PISCES [On]
 MMSPICE2 [On]
 O MMSPICE1 [On]
 PISCES [Off]
 MMSPICE2 [Off]
 j MMSPICE1 [Off]
G
 8
p^ /
. . . . . . . . .
MMSPICE1 simulations which are also shown. Some discrepancy
in the submicron region is apparent. This could be due to a
parasitic peripheralregion transistor unaccounted for in
MMSPICE2 simulations; or possibly to slightly different
physical model parameters, e.g., mobility, assumed by PISCES
and MMSPICE2. In the switchoff case, the crowding is seen
to be insignificant as implied previously. It can be
inferred then that the reduction of switchoff delay of a
scaled device is primarily caused by the reduced charge
storage rather than the diminished crowding in the BJT.
Additional verification of the NQS crowding formalism in
MMSPICE2 is provided in Fig. 3.11 where switchon transient
collector currents predicted by PISCES, MMSPICE2, and
MMSPICE1 are plotted. These currents were taken from the
WE=l 4.m simulations of Fig. 3.10. Note the good
correspondence in time between the PISCES and MMSPICE2
currents, which are separated from the MMSPICE1 current by a
significant (added NQS) delay.
In MMSPICE1, a semiempirical accounting for current
crowding can be effected by using a parameter which defines
the intrinsic base resistance as a function of the current
dependent charge. Although the parameter could account for
the current crowding for given device dimension, it is not
applicable to other device dimensions since the parameter is
neither scalable nor predictable. Hence it cannot yield a
102
13
107
0 11010 21010 31010
Time [sec]
Fig. 3.11
Predicted switchon transient collector currents
taken from the PISCES, MMSPICE2, and MMSPICE1
simulations of Fig. 3.10 for WE=1.4(Jm.
41010
trend like Fig. 3.10.
Finally, to emphasize the mixedmode NQS simulation
capability of MMSPICE2, transient simulations of an ECL
inverter stage, the basic building block of highspeed
digital circuits, were done. Fig. 3.12 shows the circuit
diagram; the four nominal BJTs have WE=1.2pm. The output
voltage waveforms of the circuit predicted with and without
(via MMSPICE1) current crowding are plotted in Fig. 3.13.
The effect of the NQS current crowding is apparent; the
propagation delay is increased by almost 50%.
3.5 Summary
A novel NQS model for transient current crowding in
advanced BJTs has been developed. The new model, based on
the use of the previous timestep solution in the current
timestep analysis, characterizes a timedependent effective
bias on the emitterbase junction for each circuit nodal
iteration at each timestep in a seminumerical analysis
following Hauser [Hau64], but physically accounting for base
conductivity modulation and the NQS nature of the crowding.
The NQS model, implemented in MMSPICE2, enables a semi
numerical, scalable, mixedmode device/circuit simulation
capability for applicationspecific TCAD. The tool is
supported by numerical simulations of advanced BJT structures
Vin O
2.5 V
Fig. 3.12 A
f
GND
Vout
5.2 V
n advancedtechnology ECL inverter circuit. The
our BJTs have LE/WE=9.22lm/l.2j2m.
0.5
0.7
0.9
1.1
1.3
1.5
0 1 1010 21010 31010 41010
Time [sec]
Fig. 3.13
Switching waveforms of the ECL inverter circuit
simulated with (MMSPICE2) and without (MMSPICE1)
the transient current crowding accounted for.
51010
using PISCES. From the simulations of a representative BJT
inverter circuit, the following conclusions were reached.
(1) For the switchon transient, peripheralemitter crowding
causes an added delay, and tends to become insignificant
only when WE is scaled to deepsubmicron values.
(2) For the switchoff transient, the added delay due to
centralemitter crowding is negligible, at least for
WE<2pm. Indeed the reduced delay of a scaled (WE and AE)
device is due predominantly to the reduced charge storage
in the BJT.
We note that the novel modeling/implementation involving
use of the previous timestep solution to update the model
for the current timestep analysis could be a viable means of
accounting for general NQS behavior in seminumerical
transient device/circuit simulation. Such behavior must
indeed be modeled to enable truly predictive mixedmode
simulation for TCAD.
CHAPTER 4
ANALYTIC ACCOUNTING FOR CARRIER VELOCITY OVERSHOOT
4.1 Introduction
In advanced siliconbased bipolar technology, the
vertical as well as the lateral dimensions of the BJT are
being scaled to deepsubmicron values. Consequently, very
high electric fields and field gradients are not uncommon in
the scaled device. When the field increases rapidly over
distances comparable to the energyrelaxation mean free path,
carrier velocity can overshoot the value corresponding to the
local electric field. This enhanced transport occurs because
the carrier (kinetic) energy, which controls the collision
time and hence limits the velocity, lags the field and
remains relatively small [Ruc72]. Such a nonlocal effect
has been recognized as significant in MOSFETs and MESFETs for
years, but only now is its significance in advanced bipolar
transistors (BJTs) becoming an issue [Lee89, Cra90].
Recent work [Fus92] has indicated that velocity
overshoot in scaled silicon BJTs can be beneficial, and must
be accounted for in the device and circuit design. The
effect, however, has not yet been physically accounted for in
any circuit simulator. Indeed, this phenomenon is not
accounted for in most device simulators because of the
implied computational intensiveness. The conventional drift
diffusion current equation used in ordinary circuit and
device simulators does not account for the nonlocal effect
of an inhomogeneous electric field on the carrier velocity.
It is based on the assumption that the drift velocity is a
function of the local electric field, and ignores the actual
dependence (of mobility) on carrier energy.
Nonlocal effects on carrier transport have been
accounted for using different analyses, but with severe
restrictions because of the accuracy/computational efficiency
tradeoff. Hence these analyseswhich include rigorous
Monte Carlo statistical treatments [Lee89], less complex
solutions of the hydrodynamic equations involving the
solution of the moments of Boltzmann transport equation
(i.e., a set of equations describing conservation of particle
number, momentum, and energy solved in conjunction with
Poisson's equation) [Blo70], and even simpler solutions of
the energy transport equations which, with some assumptions,
can be derived from the hydrodynamic model [Bor91]have
limited utility for device simulation and virtually no use
for circuit simulation. Alternatively, the socalled
augmented driftdiffusion (ADD) transport model [Tho82],
which retains most of the efficiency of the driftdiffusion
equation but uses additional analytic terms to account for
the nonlocal effects, has been proposed as a way of
efficiently extending the utility of drift/diffusionbased
tools for scaled technologies.
In Section 4.2, a simple but physical analytic model for
firstorder accounting of the electron velocity overshoot in
advanced siliconbased BJT "circuit simulation" is presented.
The model, which characterizes the nonlocal electron
velocity in the highfield collector spacecharge regions
(SCRs), is shown to be identical to the ADD formalism when
the electron diffusion is negligible. The associated
velocity relaxation, which is not accounted for in the ADD
model, is characterized phenomenologically to be consistent
with the overshoot analysis. In Section 4.3, the comparison
of our model with the energy transport analysis is presented.
In Section 4.4, the implementation of the model in MMSPICE is
discussed. In the last section, device and circuit
simulation results are presented to assess the significance
of the electron velocity overshoot in advanced silicon
bipolar and BiCMOS technologies, and to support the model.
This is the first time that a nonlocal effect has been
explicitly accounted for in a circuit simulator.
4.2 Model Development
4.2.1 Velocity Overshoot
When the randomly moving conductionband electrons in a
semiconductor encounter an electric field, they experience an
increase in average (drift) velocity, and an increase in
average kinetic energy which however tends to lag the drift
velocity [Ruc72]. When the kinetic energy is important
(i.e., when the electrons are not in thermal balance with the
lattice), a phenomenological force acting on the electrons
can be expressed in one dimension as
qE(eff) d (EC W (4.1)
dx
where Ec and W are the (average) potential and kinetic
energies of the electrons respectively. Note that Ec and W
in (4.1) are "correlated" in accord with electron flow. When
Wt is small (=3kT/2 where T is the lattice temperature),
E(eff) is the actual field, E, proportional to dEc/dx as it is
classically expressed.
Ballistic transport of the electrons, driven by E, would
result in unlimited 7W. However the electrons in a crystal
lattice frequently collide with impurities and phonons, the
result of which is to randomize their motion and limit their
(average) drift velocity, v, and hence their momentum.
Effectively the collisions give rise to a retarding force
proportional to the velocity, as characterized by the balance
of momentum [Shu81]:
m*d = qE(eff) m* v (4.2)
dt T(71
where m* is the effective mass of conduction (sub)band
electrons and T(7'4 is an energydependent momentum relaxation
time. Combining (4.1) and (4.2) yields
*dv dEc dW v
m +  m
dt dx dx (74
= qE + m (4.3)
For dc or quasistatic analysis, dv/dt=0 in (4.3) and
v m qE + (4.4)
Note that when dW/dx is negligible, (4.4) becomes a well
known equation defining the electron mobility 1(7'4 (=1v/El):
( q (4O (4.5)
m
The mobility is expressed as a function of 4' to emphasize
that it depends more on the local carrier energy than on the
local electric field. Using (4.5) in (4.4) with the chain
rule for differentiation gives
LIq dlEl E dxH
= vo(E) + L(E) dE (4.6)
L E dxJ (4.6)
where vo(E) is the conventional drift velocity defined by the
local field, and L(E)E(d'W/dlE )/q is a phenomenological
length coefficient [Pri88], which describes to firstorder
the nonlocal effect of the electric field gradient on v.
For L(E)#0, a large dE/dx in (4.6) implies a possibly
significant velocity overshoot, Ivi>vo(E) in accord with
the more rigorous physics underlying the electron transport.
Note that (4.6) is identical with the ADD formalism [Tho82]
when the diffusion of carriers is negligible [Kan91]. The
field gradient in (4.6) was substituted with the quasiFermi
level by other authors [Kiz89], to avoid inappropriate
overshoot corrections in the presence of builtin electric
field. However this would not be important in real
applications, since the simulation of the equilibrium
condition is not needed in most cases.
The length coefficient has been characterized via Monte
Carlo analysis [Art88] by several investigators. However the
results show some quantitative differences, possibly because
of the different transport parameters and band structures
used. Recently, Chen et al [Che91] derived an analytic
formula for L(E), but its utility is subject to uncertainties
in the evaluation of some model parameters. Hence we suggest
a simplified piecewiselinear representation of L(E) for
electrons in silicon at room temperature, based on Artaki's
Monte Carlo simulations [Art88], which is illustrated in Fig.
4.1. In fact, L(E) can be negative for low E, although the
velocity undershoot thereby implied by (4.6) is generally not
significant [Lun90] and will be neglected here.
Equation (4.5) implies that the classical mobility
decreases with increasing electric field since the electrons
gain kinetic energy which reduces the average (scattering)
time between collisions. When the velocity imparted to an
electron by the applied field is much less than the random
thermal velocity, T is however insensitive to E, implying a
linear v(E) dependence: vo=p1oE where 4o is the lowfield
mobility. At high fields however, the drift velocity becomes
comparable to the random thermal velocity, and T is reduced.
The drift velocity (magnitude) in this case, in the absence
of a high gradient of E, approaches a limiting (saturated)
1 105 I I I I
8106
6106
0 O
4106O O
2106
0 O
0 O Artaki's Work
2106 Our Model
4106
0 20 40 60 80 100 120 140
Electric Field Magnitude [KV/cm]
Fig. 4.1 The length coefficient versus electric field
(magnitude) for silicon at room temperature. The
points were derived from Monte Carlo simulations
[Art88], and the piecewiselinear approximation is
used in our model.
value vs (=107cm/sec in silicon at room temperature), which
can be empirically expressed as the product of o0 and a
critical electric field (magnitude) Es defining the onset of
velocity saturation: vs=oEs.
Hence depending on the magnitude of the electric field
in a region with dlEl/dx > 0, the magnitude of the carrier
drift velocity in (4.6) can be expressed as
v, = olEI [1 + L(E) dE] = olEI for IE < Es (4.7)
L E dxJ
and
= vs + L(E) dE] for IE > Es. (4.8)
L E dxJ
The typical value of Es for electrons in silicon at room
temperature is less than 30KV/cm, and for IEI
vanishes as shown in Fig. 4.1. Hence as indicated in (4.7),
Ivl=ol IEl for this case, in accord with the conventional
characterization. This simplification means that the
velocity overshoot characterization is needed only when
IEI>Es as in (4.8), and that otherwise the conventional
driftdiffusion formalism with (4.7) is still applicable even
though dlEl/dx is high.
4.2.2 Velocity Relaxation
The analytic velocity overshoot characterization in
(4.8) is strictly valid only when the magnitude of the
electric field is increasing in the drift current direction.
It would yield no overshoot when diE/dx = 0 or an undershoot
when dlEl/dx < 0, independent of the history of the
transport, and hence is nonphysical for these cases. For
example, a hot (high'4t electron entering such a region where
dlEl/dx is not positive must travel a few mean free paths to
reach the velocity corresponding to the local field, and
hence would experience velocity overshoot. This relaxation
can be neglected for MOSFETs and MESFETs because the only
significant nonlocal effects occur under the gate where
electrons are accelerated to the drain by a high field with
dlEl/dx > 0 [Kiz89, Kan91]. However for the BJT, which
contains significant (spacecharge) regions with diE/dx < 0
adjacent to those with dlEi/dx > 0, the velocity relaxation
following overshoot must be simulated. Details on various
types of SCRs will be presented in next section.
To understand the velocity relaxation in the advanced
BJT, consider a mental experiment. Fig. 4.2 shows the
possible relaxation of the drift velocity in the collector
side of the basecollector junction SCR where ]El is
v(x)
v(0)
v(O) Case 1
Case 3
vs .
Case 2
Edge of SCR
Fig. 4.2 Possible distributions of the drift velocity when
IEl is decreasing with distance. Note that the
electric field magnitude at the edge of SCR is
assumed to be Es.
decreasing with distance (see Fig. 4.3(a)). Note that the
electric field magnitude at the (nebulous) edge of the SCR is
implicitly assumed to be Es [Jeo89] Normally when a hot
electron leaves a highfield region, its velocity will
decrease with distance due to the scattering by which it
transfers its energy to the lattice (see Case 1 in Fig. 4.2).
The relaxation however becomes somewhat different when the
width of the SCR gets smaller. At a glance, it seems likely
that the velocity would not decrease very much from its value
at the junction because of the reduced scattering. But
actually this tendency would be compensated by the velocity
undershoot tendency [Lun90], which obtains when the electric
field is decreasing very rapidly. The kinetic energy
responds to fields more slowly than does the carrier
velocity; hence immediately after the high to lowfield
transition, the carrier's kinetic energy is still high, and
thus its mobility is lower than that corresponding to thermal
balance between the carrier and the lattice. After the
electron has dissipated its excess energy, it would then have
the velocity vs (see Case 2). This is supported by the fact
that L(E) in (4.8) is 0 at the edge of the SCR because El is
assumed to be Es. Taking these two conflicting phenomena
into consideration, we assume that the velocity would decay
monotonically with distance and finally reach vs at the edge
of the SCR (as described by Case 3).
Based on this insight, we use a phenomenological
representation of the velocity relaxation in an SCR where
dlEl/dx < 0 by simplifying (4.2) to
dv dv v
dt dx t
or
dv 1 ,_ v (4.9)
dx T s
where s is an average mean free path for velocity relaxation.
The solution of (4.9) is
v(x) = v(0) exp(x/s) (4.10)
where v(0) is the velocity at the point where IEl is maximum
in the SCR. Since the velocity must be continuous, v(0) is
derived from the analysis of the velocity (overshoot) in the
adjacent region where dlEl/dx > 0. To estimate s, we assume
as discussed above that the carrier velocity reaches vs at
the edge of the SCR. (This assumption is consistent with a
common designation of an SCR [Jeo89].) Thus
s WRR (4.11)
InV(O)
Vs,
where WRR is the width of the relaxation region.
4.2.3 Effective Saturated Drift Velocity
To this point, we have modeled the hotelectron velocity
in an SCR using either the length coefficient or the
scattering mean free path, depending on the sign of diE/dx.
To facilitate the implementation (discussed later) of the
model into the bipolar device/circuit simulator MMSPICE, we
define now an effective saturated drift velocity vs(eff) based
on the actual transit time of electrons in the SCR being
analyzed:
dt = WSCR (4.12)
IC I C v(x) Vs(eff)
WSCR WsCR
where v(x) is given by (4.8) or (4.10), and WSCR is the width
of the SCR in which IEl is greater than Es.
For the advanced BJT, different operating conditions are
distinguished by the charge conditions [Jeo89] in the
epitaxial collector region, as reflected in Fig. 4.3. The
electric field distributions shown are determined by the bias
on the basecollector junction and the collector current.
Wvs 0 Wscc
E(x)
0 WQNR
WEPI WBL
E(x)
WEPI
E(x)
Fig. 4.3 Electric field distributions in a basecollector
junction SCR (a), and a currentinduced SCR (b)
associated with nonohmic quasisaturation, i.e.,
base pushout. When either SCR expands, the entire
epi layer can become spacecharged (c).
Fig. 4.3(a) represents the conventional junction SCR at the
basecollector junction under lowcurrent conditions. For
highcurrent conditions, when nonohmic quasisaturation
(base pushout) prevails, the currentinduced SCR exists in
the epicollector as denoted in Fig. 4.3(b). Note that the
electric field is assumed to be Es at the edge of the
collectorside SCR in both cases; this assumption in fact
defines the SCRs [Jeo89] When either SCR expands, the
entire epi layer can become spacecharged, as shown in Fig.
4.3(c). We must consider the three SCR types in the BJT
separately.
4.2.3.1 Junction SCR
When the SCR exists across the basecollector junction,
as shown in Fig. 4.3(a), (4.12) applied to it yields
f0 fWscc
dx + dx Wvs + Wscc (4.13)
V(x) v(x) Vs(eff)
where Wvs and Wscc are the widths of the base and collector
sides of the SCR respectively. The carrier velocity v(x) is
evaluated depending on the sign of dIE/dx. In the base
side, the velocity is characterized via (4.8), using the
depletion approximation coupled to a firstpass (vs(eff)>Vs)
MMSPICE simulation to describe E(x) and Wvs:
dE q [NA(x) + n]
dx E
=  NA(X)
 NAO exp( x + WBM) (4.14)
C LWBM
where the assumed exponential doping profile is consistent
with the basetransport analysis of the BJT model [Jeo89] in
MMSPICE; WBM is the metallurgical base width. Thus
(4.15)
E(x) qWBM NAO exp +x +W C (4.15)
EM [WBM
The integration constant C can be easily evaluated from the
electric field at the junction (x=0), which is available from
the output of the BJT model routine in MMSPICE. E(x) and
dE/dx are then substituted into (4.8) to give v(x) for the
first integral in (4.13). The validity of using the
depletion approximation here will be discussed in Appendix D.
In the collector side, (4.10) is used directly for the
second integral in (4.13), with v(0) being equated to that
derived from the analysis of base side. Both integrals in
(4.13) are evaluated by a numerical method to give vs(eff).
Strictly, the value of Es in the base side tends to be
greater than that in the collector side because the electron
mobility (Io) in the base is lower due to the higher doping
concentration. However because the (compensated) doping is
generally not known precisely and because this variation in
Es is only a secondorder effect, we neglect it.
4.2.3.2 Currentinduced SCR
When the currentinduced SCR exists, as illustrated in
Fig. 4.3(b), (4.12) applied to it yields
WEP dx WEPI WQNR (4.16)
d ,W, (4.16)
v(x) Vs(eff)
where WQNR is the extended width of the pushedout (quasi
neutral) base region. The transit time across the portion of
the SCR in the adjacent buried layer of the BJT structure is
neglected since the heavy doping there implies only a
negligibly thin depletionregion width, WBLWEPI.
From the firstpass MMSPICE simulation [Jeo89], the
electric field in the SCR and WQNR are obtained in accord with
dE An
dx E
C \qAvs I
where NEPI is the doping concentration of the epicollector
layer and An is the excess electron density in the SCR, which
is assumed to be spatially constant since the current Ic is
constant. From (4.17),
E(x) = q Ic N (x WQNR) Es (4.18)
Equations (4.17) and (4.18) are substituted into (4.8) to
yield v(x), and vs(eff) is evaluated from (4.16). In the
vicinity of the boundary between the SCR and the quasi
neutral region in the epicollector, the electric field
gradient is very large. However this transitional region can
be ignored because the length coefficient is, as shown in
Fig. 4.1, assumed to be 0 when IEI
4.2.3.3 Special case
In previous sections, the velocity overshoot effect was
characterized via regional analyses depending on the sign of
the field gradient. There is a special case for the BJT
however where the overshoot effect would not be properly
accounted for in this manner. This is the case where the epi
layer is completely spacecharged, and the magnitude of the
electric field is still increasing with distance due to non
ohmic quasisaturation, as shown in Fig. 4.3(c). (Note that
when the entire epi layer is spacecharged, but El is
decreasing with distance, the overshoot analysis for the
junction SCR is still applicable.)
According to our formalism, the same overshoot analysis
would be applied in the collector side as in the base side.
Of course, this is adequate if the field gradient is
relatively large. When the electric field is increasing
slightly however, the direct application of our model would
tend to exaggerate the overshoot effect since velocity
relaxation is ignored. In fact, the carrier velocity would
decrease with distance in the epicollector. In order to
cope with this deficiency of our formalism, we empirically
combine the overshoot model with the relaxation model for
this case as follows:
v(X) = Vrei(X) + ov (x) exp  + Voffset (4.19)
f dx)
where vrel(x) and vov(x) are the velocity distributions
characterized by the relaxation and overshoot models
respectively, and f is an empirical weighting factor. When
the field gradient is very small, (4.19) reduces to (4.10),
implying that velocity relaxation would be predominant in the
collector side. When the gradient becomes large, v(x) is
given as the sum of vrel(x) and vov(x) with the empirical
factor chosen to ensure a smooth transition from velocity
relaxation to velocity overshoot. The offset velocity,
offset in (4.19) is used to make the velocity at the junction
continuous.
4.3 Comparisons with Energy Transport Model
One way to characterize the velocity overshoot effect is
to solve the energy transport equation [Bor91]. Such a
solution can provide support for our simple analytic model.
In this section, we will numerically solve the energy
transport equation Goldsman et al presented [Gol88], and
contrast it with our model. By assuming the electron energy
as entirely thermal, they derived the steadystate momentum
equation from the Boltzmann transport equation as
S= (w) qE 2 dw 2w dn (4.20)
m* 3 dx 3n dxx
1 3 3
where w is the average electron energy (= m*v2 + kTe = kTe
where Te is the electron temperature), Tp(w) is the energy
dependent momentum relaxation time, and n is the electron
concentration. Combining (4.20) (with dn/dx=0) with the
steadystate energy equation, Goldsman et al derived an
equation for average electron energy that includes the effect
of velocity overshoot:
dw 21 qE 9 40 m* (o) + 1/2 (4.21)
dx 20 20 [9 Tpw qE2]
where w (w) is the energy relaxation time, and wo is the
thermal energy of the lattice (=3kT/2).
In order to solve these equations, both Tp and Tw must be
known as functions of the electron energy. Although Goldsman
et al evaluated the relaxation times by Monte Carlo
simulations in homogeneous fields, we use simple functions to
empirically approximate the parameters they derived:
Tp(w) = co + c and (4.22)
w
Tw(w) = do + d1w + d2w2 + d3w3 (4.23)
where Cn and dn denote empirical constants. In Fig. 4.4, the
discrete points represent the momentum and energy relaxation
times Goldsman et al have derived, and the solid lines which
best fit the data are given by (4.22) and (4.23). Then the
energy dependent carrier velocity can be numerically
evaluated from (4.20) and (4.21), since those equations are a
function of the single variable w.
0.08
0.06
0.04
0.02
0
0 0.1
Fig. 4.4
0.3
m
r3
CD
0.28 <
CD
o
I
0.26 3
CD
"0
Zcn
0.24
0.2 0.3 0.4 0.5 0.6 0.7 0.8
Average Electron Energy [eV]
Momentum and energy relaxation times as functions
of energy.
For comparisons, we evaluated the velocity distributions
for the typical advanced BJT, when the junction or the
currentinduced SCR exists, using our model and that of
Goldsman et al. Fig. 4.5(a) shows the predicted velocity
distributions in the junction SCR when VBE=0.7V and VBC=O.OV
are applied to the terminals of the device. (For the
effective mass of conduction subband electron, m*=0.26mo was
used, where mo is the rest mass [Mu189]. The electric field
used as inputs for both the models was available from the
output of MMSPICE.) As described before, our overshoot
analysis is done when the magnitude of the electric field is
increasing (x<0). In accord with our piecewiselinear L(E)
model, the carrier velocity reaches its peak value Vpeak when
the length coefficient is at its maximum value at E=50KV/cm
(see Fig. 4.1). Note that the location of peak is about the
same as that predicted by the energy transport model. When
reverse bias is applied on the basecollector junction (VBC=
2.0V), peak increases as shown in Fig. 4.5(b), because the
gradient of the electric field also increases. Figs. 4.6(a)
and (b) illustrate the velocity distributions in the current
induced SCR (for VBE=1.OV, VBC=O.OV or VBC=2.0V) .
We note in the above figures that our model predicts a
higher peak overshoot velocity than that yielded by either
the energy transport model or Monte Carlo simulations (not
88
(a)
3.5 107 i I ,I I I 140
7
3.0 10 Our model 120
S  Energy transport
7
2.5 10 100
7
Q 2.0 10 80
. 
7 <
1. 510A A60 
1.0 10  40
5.0 106 20
0.1 0 0.1 0.2 0.3 0.4
Distance from BC junction [um]
(b)
3.5107 I I I i 140
S3.0 Our model
3.0 10 120
S   Energy transport
7 !
2.5 10 100
N 10
0) 7
U 2.010 80 x
7 <
S1.510 6
1.0 107   40
I  \_
5.0106 20
0.1 0 0.1 0.2 0.3 0.4
Distance from BC junction [urn]
Fig. 4.5 Drift velocity and electric field in junction SCR:
(a) VBE=0.7V and VBC=O.OV; (b) VBE=0.7V and VBC=
2.0V.
2.0 107
1.5 107
1.0 107
5.0 106
0
2.0 107
1.5 107
1.0 107
5.0 106
, I I I I ,I
0 0.1 0.2 0.3 0.4 0.5
Distance from BC junction [urn]
(b)
.I.I I I .1 1 1 1 1 1 1 1 1 I ,
l r I I I I I I I I I I I~l 1 I I I I r
 200
 150
 100 0
50
250
 200
m
 150 "
 100 3
50
Distance from BC junction [um]
Fig. 4.6 Drift velocity and electric field in current
induced SCR: (a) VBE=I.OV and VBC=0.OV; (b) VBE=1.0V
and VBC=2.0V.
 Our model
   Energy transport
/"
 Our model
   Energy transport
I
I" /
/
/
/
I
/
/
. . I I, , I
shown in the figures) [Prof. M. Lundstrom of Purdue
University, private communication, 1991]. This discrepancy
could mean that the length coefficient [Art88] which we used
might be erroneous. Indeed when the highIEI saturation
value of L(E) in Fig. 4.1 is reduced from 4.5x106cm to
2.0x106cm, which has been suggested [Art88], the MMSPICE
predicted velocity overshoot is in better agreement with that
predicted by the energy transport and Monte Carlo analyses.
This uncertainty in L(E) can be attributed to the different
set of transport parameters used. However we stress that the
terminal characteristics of advanced BJTs predicted by our
model, which will be shown later, agree quite well with
results [Fus92] of measurements and numerical simulations
based on a hydrodynamic model for energy transport.
Conversely then, we note that the energy transport model has
several uncertainties as well. It is based on several
equivocal assumptions. For example, it assumes that the
electron energy is entirely thermal. Also, the results
depend on the degree of the energy transport equation, and
there are still some uncertainties in the evaluation of the
model parameters such as m*, Tp, and T,. Monte Carlo analysis
is not unequivocal either. For example, detailed and
accurate information about the numerous scattering parameters
as well as needed details of the energyband structure are
lacking. With these deficiencies then, our model is
reasonable for firstorder accounting of the electron
velocity overshoot in circuit simulation, which has never
been done before.
4.4 Implementation
The implementation of the electron velocity overshoot
model in MMSPICE is based on a single iteration of the
existing (conventional) model routine [Jeo89] for the (n+pnn+)
BJT, as illustrated in Fig. 4.7. The analysis is done for
each iteration of the circuit nodal analysis at each time
step. With VBE and VBC passed in from the nodal analysis, the
(onedimensional) BJT model routine, which assumes a
saturated drift velocity vs (no overshoot) in the highE epi
collector SCRs, is called to solve the conventional ambipolar
transport, and characterize E(x) Thus unlike empirical
circuit models, the MMSPICE BJT model is susceptible to an
extension to account for the augmented nonlocal carrier
velocity distribution. From the predicted E(x), combined
with the length coefficient L(E), the carrier velocity is
evaluated depending on the SCR type (see Fig. 4.3). Then
from v[E(x)], the effective saturated drift velocity,
vs(eff)>Vs, is evaluated as described in Section 4.2.3.
Fig. 4.7 Algorithm for implementation of velocity overshoot
model in MMSPICE.
Once Vs(eff) is characterized, Es is correspondingly
updated as well to Es(eff)yVs(eff)/Io, which is higher than the
preliminary value. Fig. 4.8 illustrates the resulting
velocityfield model in the epicollector SCRs used in
MMSPICE, which we believe is suitable for firstorder
accounting of the velocity overshoot in circuit simulation.
Note that the v(E) slope (i.e., the lowfield mobility p0) is
not changed. Hence the solution obtained in regions where
IEI
physically appropriate.
With vs(eff) and Es(eff), the MMSPICE BJT model routine is
called once more to effect the firstorder accounting for the
nonlocal transport in the predicted device currents and
charges. The accounting for velocity overshoot, which is
done here in a circuit simulator for the first time, is
computationally efficient, and enables representative mixed
mode simulation for advanced bipolar technologies as we now
demonstrate.
4.5 Simulations
In this section, MMSPICE device and circuit simulation
results are presented to assess significance of the velocity
overshoot effects in advanced siliconbased bipolar
