MicroelectronicEngineering19 (1992) 31-38 Elsevier
Simulation
of Compound
Semiconductor
31
Devices
Wim Schoenmaker and Rudi Vankemmel IMEC Kapeldreef 75, B-3001 Leuven-Heverlee, Belgium Abstract We r e v i e w t h e p r e s e n t s t a t u s of compound s e m i c o n d u c t o r device simulations. Some scaling techniques are discussed for the simulation of abrupt heterojunctions. An example illustrates the method. 1. INTRODUCTION Although Silicon technology represents still the bulk of the semiconductor industry, compound semiconductor devices have a t t r a c t e d m u c h interest. There exists several reasons for it. First of all, processing facilities, have made it possible to fabricate devices with a much higher precision as a few years ago. Secondly there is an increasing demand for devices with very high cutoff frequencies a nd fast switching times. Therefore in m a n y laboratories, engineers are experimenting on heterostructure devices which m a y consist of m a n y layers of materials. It would be of great help, if one could perform reliable computer simulations of these devices in order to speed up these experiments and moreover in order to understand better the interplay between the electronic behaviour and the presence of the various layers with different material properties. Unfortunately, we are still far removed from the situation that we can offer a simulation tool that would satisfy the needs of the engineers in the compound-semiconductor industry, in comparison with their colleagues from the Silicon industry. In the latter one, the designers of CAD tools have shifted the emphasis of the new developments away from the simulators towards the integration of process- and device simulators, complemented with advanced a u t o m a t i z e d optimalization algorithms, i.e. c o m p u t e r - i n t e g r a t e d manufacturing (CIM) [1]. Beside the integration, let us summarize the requirements which a device si mul a t or for compound-semiconductor devices should fulfil. These are: robustness, flexibility, reliability and user-friendliness. These demands are far more difficult to realize for compound-semiconductor device simulators, t h a n for Silicon device simulators, as will become clear in the sequel of this paper. 2. SURVEY OF RECENT RESULTS The simulation strategies that one can option for may be subdivided into three mainstreams, albeit that the splitting remains somewhat artificial, since results and methods from one approach often serve as input for a n o t h e r one. In the Monte Carlo approach, one simulates the stochastic flow of the particles through .the device by mimicing this flow process by using r a n d o m variables. 0167-9317/92/$05.00 © 1992 - ElsevierSciencePublishersB.V. All rights reserved.
32
W. Schoenmaker, R. Vankemmel / Simulation of compound semiconductor devices
In order to determine the electric force that the particles experience, Poisson's equation is solved in a deterministic way. This brings us to the second approach, which consists of the formulation of the balance equations for the carrier densities and energy densities and solve t h e s e t o g e t h e r w i t h Poisson's e q u a t i o n , all in a d e t e r m i n i s t i c way (hydrodynamical modeling). The Monte Carlo method will be discussed by P. Lugli, [2] and we will not consider it here any further. Beside the hydrodynamical approach, quasi-lD models have been developed [3] for a n u m b e r of structures. Such models are important because of their low CPU cost such t h a t they can be easily incorporated into circuit simulators. Elaborate simulations of GaAs-ALGaAs MODFET structures have been performed by Shawki, Salmer and E1-Sayed [4,5]. The simulations are based on a 2D hydrodynamical model and include hot-electron effects and the presence of DX centres. The inclusion of the DX centres in the simulation showed t h a t they strongly reduce the cutoff frequency. The heterojunction of the GaAsA1GaAs interface was included by grading the conduction band discontinuity of .23 eV over 20 .~. A t h r e e dimensional device simulator for GaAs-A1GaAs heterojunction devices was developed by Chan and Shieh [6]. The simulator essentially solves the drift-diffusion moflel in steady state and the heterojunctions are graded from 100 A to 1000 A. Since the minimal grid spacing is 100 A, the 100 Agrading is considered to be an abrupt junction. A n o t h e r s i m u l a t o r tool for GaAs-A1GaAs h e t e r o s t r u c t u r e s has been developed by HITACHI. It is an extended 2D device simulator based on the drift-diffusion model in which a single donor level DX centre model was i n c l u d e d [7,8]. The s i m u l a t i o n s confirm th e l a r g e d e c r e a s e of the transconductance due to the existence of the DX centres. For a n u m b e r of questions concerning heterojunction devices, such as the c a r r i e r d i s t r i b u t i o n at the hetero-interface one m u s t a d d r e s s q u a n t u m mechanical methods. This implies that SchrSdinger's equation and Poisson's equation have to be solved self-consistently. This problem has been investigated by several groups. Recently, progress has been made also, to describe t r a n s p o r t along t he s e lines. Magnus et al. [9] choose a l i m i t e d set of parametrized base functions in order to describe the subbands, in which the p a r a m e t e r s have to be determined self-consistently with the solution of the Poisson equation. The qua nt um transport problem, lateral to the interface, can be solved by using the m e a n field approximation and the quantum-balance equations for m o m e n t u m and energy [10,11]. The same problem has also recently addressed by Takano, Yu and Dutton [12]. These authors approach the quantization problem by dividing the device in q u a n t u m - m e c h a n i c a l and classical regions. At the interface, q u a n t u m effects dominate the carrier distribution and a window selects the region where SchrSdinger's equation is solved. Such an approach limits the use of quantummechanical methods to regions where their outcome differs from classical methods, but a fully numerical approach requires a careful t r e a t m e n t of the boundary conditions. So far, one has been working with plane waves (nonconfined electrons) or exponentially decaying waves (confined electrons) to match the wavefunction outside the window.
W. Schoenmaker, R. Vankemmel / Simulation of compound semiconductor devices
33
3. ABRUPT HETEROJUNCTIONS The simulation of abrupt heterojunction is usually done by grading it over a minimal mesh distance. Such an approach hampers the use of the simulator by external users because one may easily take the grading too steep such that one pushes into an instability. MBE techniques allow the deposition of single atomic layers (-6 /k) which would require extremely fine m e s h e s at the interface. In order to avoid the problem due to too-steep grading, we have idealised the heterojunction by a two-valued entity. In Fig. 1 a schematic view is given of a heterojunction. Since the electrostatic potential is by definition sited at the middle of the bandgap it becomes a two-valued function at the heterojunction nodes. Only one value is stored in the computer memory whereas the other one is calculated from the bandgap discontinuity, if required. In this way we could construct a compound-semiconductor device simulator, PRISM [13,14], which does not require detailed knowledge of the interplay between mesh generation and the presence of heterojunctions (user-friendlyness). Another important issue of abrupt heterojunctions concerns the t r e a t m e n t of the current continuity equations at the junction nodes. |
J.4 I
MAT UP
~1 space charge region
I
ap_down
Ei E
J (1-f).E diff gap_up
Ec
MAT_DOWN Ediff = E mat_up- E mat_down
Figure 1. Schematic band diagram at a heterojunction. These nodes have also an abrupt jump in the carrier concentration when crossing the heterojunction interface. Although it is well possible to perform the assembling of the decoupled continuity equations, it becomes ambiguous to assign the carrier concentration update charge either to the upper or lower concentration at the junction node. However, when using a coupled scheme the ambiguity can be avoided, because the update does not occur on the charge variable but the fermilevel, which is considered to be continuous when crossing the junction. A l t h o u g h h e t e r o j u n c t i o n transistors are majority carrier
34
W. Schoenmaker, R. Vankemmel / Simulation of compound semiconductor devices
(electrons) devices, it might be important to simulate the influence of p-type GaAs buffer layers, which can only be done accurately by solving the hole continuity equation [15]. Moreover, the 2DEG-HBT [8] is a bipolar device which makes the hole equation indispensable. The simulation of heterojunction devices in bipolar mode requires a careful scaling algorithm. First of all, the conventional de Mari scaling is applied. This gets rid of unnecessary numerical factors and dimensions. It makes Poisson's equation suitable for a numerical implementation. In order to solve the current continuity equations the box integration method is used. Applying the method to the finite element discretisation scheme as shown in Fig. 2 gives the following equation for the node i [14].
N
~Jij.dij +sq.(R-G)iA i =0 j=l in which N represents the n u m b e r of neighbouring nodes of i, Jij is the S c h a r f e t t e r - G u m m e l - T a n g current from node i towards node j, posibly extended with extra components for tunneling through the heterojunction barrier and thermionic emission. Furthermore, s=+l for holes and -1 for electrons. I
k
MATERIAL UP
m
j
heretojunction
heterojuncfion
MATERIAL DOWN
Figure 2. Graphical representation of currents out of a node. In order to bring the nodal equations to comparable numerical range each nodal equation is scaled with respect to the nodal concentration. It should be noted, t h a t this scaling factor is unique for each node, therefore for the heterojunction nodes we apply a scaling with respect to the geometrical average,
W. Schoenmaker, R. Vankemmel / Simulation of compound semiconductor devices
35
•]Cup.Cdown
where Cup and Cdown are the concentrations in the upper and lower material respectively. Another important topic is the selection of a global concentration n* i (c.f. n* i in the de Marl scaling scheme ) to which all concentrations in all m a t e r i a l s are scaled. For m a n y i n t e r e s t i n g applications, the intrinsic concentration of one material can be much smaller than for another material being present in the same device. By just taking the maximum of all intrinsic c o n c e n t r a t i o n s one can push the n u m b e r s for the m i n o r i t y carrier concentrations of the materials with small n*i below the accuracy limit of the computer, which results into numerical noise and prevents convergence. In PRISM several scaling selections are built in , for example Nmat ~ log(n*i) n*i = exp Nmat i=l although it sometimes unavoidable that the user selects a trial value for n* i (robustness).
When setting the pre-defined (default) values for the v a r i o u s m a t e r i a l parameters, it turns out that the experimental data is rather limited. This is in particular the case for the carrier mobilities, the SRH recombination times, the impact-ionization coefficients and energy relaxation times. Of course, this limits the predictive power of the simulator (reliablity). On the other hand, if steady state I-V characteristics are available, one can get an estimate of the mobilities. This has been done in [13], where the device simulator has been used to estimate the electron mobility in InGaAs. Recently, in collaboration with SILVACO Data Systems, we have improved the user-friendlyness of PRISM by linking the graphical environment of TonyPlot to PRISM. Moreover, the mesh generation and device design can be done menu-driven. Fig. 3 shows a typical screen of the mesh generator for the simulation of a heterojunction device. To summarise, we have been able to construct a flexible and user-friendly simulation package for abrupt heterojunction structures. However, the predictive power is limited due to the lack of experimental data. Moreover, the robustness remains somewhat delicate. This is intrinsic to the structures under consideration because of the wide range in intrinsic concentrations that may be present in a single compound-semiconductor device. By extended the set of linear solvers it turned out to improve the robustness considerably. 4. EXAMPLE In this section we will present the results of the simulation of a SiGe heterojunction device. The lay-out of the device is given in Fig.3. In Fig. 4,5 and 6 the potential, the electron concentration and the hole concentration are depicted under the gate (x=2.0 micron) for gate voltages of-l.0 to -5.0 Volt. The source-drain voltage is -5.0 Volt. One clearly observes the presence of the abrupt heterojunction. In Fig. 7 the IDS-VGs characteristic is shown.
36
W. Schoenmaker, R. Vankemmel / Simulation of compound semiconductor devices
8q o
~=
n=1.0E15
Silicon 4 . 0 mu
SUBSTRATE
I
Fig.3 Lay-out of SiGe device. The Germanium layer has a thickness of 0.03 micron.
O;
106
. . . . .
"
''"
'! . . . . .
~ '1' . . . . . ' 1 " ' ' ' ' "
Vg=--SV
105
0
104 0
103
I
0
100
= 10
0
1 0
0
1
si
0
0.01 t 0 -3
:
L~%r-%Qf=SEIO
0
Constant
10 -4
10 -5 10
6 ~J 0.08
RVK ii
0
09
Itl
I]1
0
III
1
I1~
0
1i
Ii
,!
0
ii
12
till
0
ill
13
R~K g2
92 ~ I~
0
n
10 14
o o8
o 09
o ':,-:sa
Fig.4
rr-~obility::
0
Electron concentration (cm-3) under gate vs depth (~m). Vd=-5V, X= 2.0~m.
Fig.5
n
I
i
o
11 i
i
n
i
o C. ~
i
,i
te
I
o
II
is
JJY~l
Hole concentration (cm-3 under gate vs depth (~m). Vd=-5V, X= 2.0#m.
t4
W. Schoenmaker, R. Vankemmel Simulation of compound semiconductor devices -02
'~
-03
:'
'''"
'I . . . .
" 'i '¸ " "
' I r ''''
']'
"''
''' . . . .
''''
' I"'"
'"~''
"
~
-~
/
-04
-0
37
5
/ //
07
/
08 --0 91 !
1 1 0 O8
Y
~JJ ~o~stcpt ~,~o~iity ='u~-, .Of =SE'O q '7< 92 ,, I~ ,,~ O
09
0
1
~Pris~-,
!~ O
11
-0.:
0
,,
12
,!; 0
,,, 13
2xlO 0
14
5 o
3
E
i
v,ss Ivl
I-r" J
Fig. 6 Poisson potential (Volts) under gate vs depth (mu). Vd=-5V, X=2.0 mu.
4
Fig.5
Drain-source current(A/#m2) vs. gate voltage (Vt curve). Vd=-0.1V, X=2.0 ~tm.
5. REFERENCES 1. R. Booth et al., proceedings Internat. Workshop on VLSI Proc. and Dev. Modeling, VPAD, Oiso Kanagawa, Japan (1991) 2. P. Lugli, contribution to ESSDERC'92 conference 3. C. Snowden, proceedings of the 5-th meeting of the GaAs Simulation group, October 12-13, 1989, Villa Gualino, Torino Italy 4. T. Shawki, G. Salmer and O. E1-Sayed, IEEE Trans. CAD 9 (1990) 1150 5. T. Shawki, G. Salmer and O. E1-Sayed, IEEE Trans. Elec. Dev.37(1990)2 6. H.-C. Chan and T.-H. Shieh, IEEE Trans. Elec Dev. 38 (1991) 2427 7. H. Mizuta et al.,IEEE Trans. Elec. Dev. 36 (1989) 2307 8. P. D. Rabinzohn et al., IEEE Trans. Electr. Dev. 38 (1991) 222 9. W. Magnds, C. Sala and K. De Meyer, Phys. Stat. Sol.(b)133 (1991) K31 10. W.Magnus, C.Sala and K. De Meyer, J. Appl.Phys. 69 (1991) 7689 11. W.Magnus, C.Sala and K. De Meyer, Phys. Rev. B43 (1991)9045 12. C. Takano, Z. Yu and R.W.Dutton, IEEE Trans. CAD 9 (1990) 1217 13. W.Schoenmaker et al., COMPEL 10 No4 (1991) 631 14. R. Vankemmel et al., Solid-St. Electron. 35 No4 (1992) 571 15. R. Vankemmel et al. submitted to IEEE Trans. CAD.