2.2 Simulation Technique
2.2.1 Physical Modelling and Numerical Methods
The technology computer-aided design (TCAD) simulations, such as process and de-vice simulations, are widely used for the analysis of semiconductor dede-vices. The process simulation can generate the device geometry and doping profile according to the parame-ters of the fabrication processes. The output of process simulation is then used in the device simulation to estimate device characteristics. The drift-diffusion (DD) and hydrodynamic (HD) models play a crucial role in the development of semiconductor device simulator in the macroscopic point of view. The DD model was derived from Maxwell’s equation as well as charges’ conservation law and has been successfully applied to study device trans-port behavior, in the past decades. It assumes local isothermal conditions and is still widely employed in semiconductor device design.
Classical drift-diffusion model consists of at least three coupled partial differential equations (PDEs) for, such as electrostatic potential and electron-hole densities. When
device channel is specified, a set of the DD equations in semiconductor device simulation
where φ is the electrostatic potential and its unit is volt. n and p are classical electron and hole concentrations (cm−3). q is the elementary charge and its unit is coulomb. The net doping concentration is D(x, y, z) = ND+(x, y, z)−NA−(x, y, z). R is the net recombination rate (cm−3s−1). The carrier’s currents densities are given by
Jn = −qµnn 5 φ + qDn5 n + unkBn 5 Tn, (2.6)
and
Jp = −qµpp 5 φ + qDp5 p + upkBp 5 Tp, (2.7)
where µnand µpare the carrier mobility (cm2/V − s). The diffusion coefficients, Dnand Dp(cm2/s), satisfy the Einstein relation.
2.2 : Simulation Technique 19
The mobility model used in the device simulation, according to Mathiessen’s rule [39–
41], can be empirically expressed as:
1
µ = G
µsurf aps + G
µsurf rs + 1
µbulk, (2.8)
where G = exp(x/lcrit), x is the distance from the interface and lcrit is a fitting parameter.
The mobility consists of three parts: (1) the surface contribution due to acoustic phonon scattering
µsurf aps = B
E + C(Ni/N0)τ
E13(T /T0)K, (2.9)
where Ni = NA+ ND, T0 = 300K, E is the transverse electric field normal to the inter-face of semiconductor and insulator, B and C are parameters which based on physically derived quantities, N0 and τ are fitting parameters, T is lattice temperature, and K is the temperature dependence of the probability of surface phonon scattering; (2) the contribu-tion attributed to surface roughness scattering is
µsurf rs = ((E/Eref)Ξ
δ +E3
η )−1, (2.10)
where Ξ = A + α(n+p)N(Ni+N1ref)υυ , Eref = 1 V /cm is a reference electric field to ensure a unitless numerator in µsurf rs, Nref = 1 cm−3is a reference doping concentration to cancel the unit of the term raised to the power υ in the denominator of Ξ, δ is a constant that depends on the details of the technology, such as oxide growth conditions, N1 = 1 cm−3, A, α and η
are fitting parameters; (3) and the bulk mobility is
µbulk = µL(T−ξ
T0 ), (2.11)
where µLis the mobility due to bulk phonon scattering and ξ is a fitting parameter.
The quantum mechanical effects should be considered in the device simulation when the dimensions of the devices shrunk into nanometer scale. Various theoretical approaches have been presented to study the quantum confinement effects, such as full quantum me-chanical model (e.g. nonequilibrium Green’s function) and quantum corrections to the clas-sical drift-diffusion (DD) or hydrodynamic (HD) transport models. A set of Schr¨odinger-Poisson (SP) equations has been applied to study the quantum effect in the inversion layers, but it is a time-consuming task in the TCAD application to realistic device characterization.
Therefore, various quantum correction models, H¨ansch, modified local density approxima-tion (MLDA), effective potential (EP), density gradient (DG) and so on, have been pro-posed for classical DD or HD transport models. In this investigation, the density gradient was coupled with the DD model and solved for the quantum mechanical effects. The den-sity gradient equation can be expressed as,
J~n= −qµnn 5 φ + qDn5 n − qnµ 5 (2bn52√
√ n
n ), (2.12)
where bn = ~2/(12qm∗n) and m∗n = mk× 9.11 × 10−31kg. bn in Eq. (2.12) is the density gradient coefficient which determines the strength of the gradient effect in the electron gas.
2.2 : Simulation Technique 21
The last term in the right hand side of Eq. (2.12) is referred to as “quantum diffusion”, which makes the electron continuity equation has a fourth-order partial differential equa-tion. Therefore, such an approach is highly sensitive to noise in the local carrier density, and the methodology is highly important in cases of strong quantization. To calculate the numerical solution of the multidimensional density-gradient model, firstly we decouple the coupled partial differential equations (PDEs); approximated with the finite volume method over nonuniform mesh. The corresponding system of the nonlinear algebraic equations is then solved with the mixed monotone iteration and Newton’s iteration methods. Iteration will be terminated and postprocesses will be performed when the specified stopping criteria for inner and outer iteration loops are satisfied, respectively.
The nominal channel doping concentrations are 1.48 × 1018cm−3 and the Vth are cal-ibrated for 16-nm-gate MOSFETs. For RDF, to consider the random fluctuation effect of the number and location of discrete channel dopants, 758 dopants are randomly generated in a large cube (80 nm × 80 nm × 80 nm), in which the equivalent doping concentration is 1.48 × 1018 cm−3, as shown in Fig. 2.4(a). To determine the location of each dopant, first we use Mersenne Twister algorithm [42] to generate double-precision values in the closed interval [2−53, 1 − 2−53]. Then multiplied the generated values by the side length of the large cube, we obtain x-, y-, and z-coordinate value of the dopant. The large cube is then partitioned into 125 sub-cubes of (16 nm × 16 nm × 16 nm). The number of
dopants may vary from zero to 14, and the average number is 6, as shown in Figs. 2.4(b) and 2.4(c), respectively. In principle, 3D device simulation with the 125 channel structures almost covers cases, shown in Fig. 2.5, and thus will be fairly meaningful to reflect statis-tical randomness of dopant number. We have noticed that in this simulation only dopant within the channel region is treated discretely.The doping concentrations remain contin-uous in the source/drain region because the volume of source/drain region is two-order magnitude greater than that of channel region. Similarly, we can obtain the distribution of dopant number for the 65-nm-gate transistor, in which the dopant number may vary from 70 to 130 as shown in Fig. 2.6. These sub-cubes are equivalently mapped into the de-vice channel for the 3D “atomistic” dede-vice simulation with discrete dopants, as shown in Fig. 2.7(a). In “atomistic” device simulation, the resolution of individual charges within a conventional drift-diffusion simulation using a fine mesh creates problems associated with singularities in the Coulomb potential [43–45]. The potential becomes too steep with fine mesh, and therefore, the majority carriers are unphysically trapped by ionized impurities, and the mobile carrier density is reduced [43–45]. Thus, the density-gradient approxima-tion is used to handle discrete charges by properly introducing related quantum-mechanical effects, and coupled with Poisson equation as well as electron-hole current continuity equa-tions [46–54]. Fig. 2.7(b) shows the studied SOI FinFET with aspect ratio (defined by the fin height / the fin width) equal to two. Without losing generality, the SOI FinFET is with
2.2 : Simulation Technique 23
16-nm-gate and 1.48 × 1018cm−3 equivalent channel doping concentration. The explored 6T and 8T SRAM circuits are illustrated in Fig. 2.8(a) and 2.8(b), respectively. As the di-mension of devices continuously scaling, the 8-transistors structure is proposed to increase the static noise margin. Unlike the 6T SRAM, the access transistors are turned off during the read operation in the 8T cell, the stored data is read out through additional transistors M1 and M2 passively, which avoid the bit line impact the data directly. Thus increase the SNM. All cell ratios (CR; CR = (W/L)driver transistor
(W/L)access transistor) of the SRAM cells in this thesis are first set at unitary. The applied voltages of 16 nm and 65 nm devices are 1.0 and 1.2 volt, re-spectively. The physical model and accuracy of such large-scale simulation approach have been quantitatively calibrated by experimentally measured results [5–9]. Similarly, we can generate 125 discrete-dopant-fluctuated cases for PMOSFET through the flow of Figs. 2.4-2.6. Then, 125 pairs of NMOSFETs and PMOSFETs are randomly selected and are used for the examination of circuit characteristics fluctuations. Furthermore, we apply the sta-tistical approach to evaluate the effect of PVE, in which the magnitude of the gate length deviation and the line edge roughness follows the projections of the ITRS 2007 [55]. The 3σ for process variation induced gate length deviation and line edge roughness are 1.5 nm and 4.3 nm for the 16 nm and 65 nm devices, respectively. In estimating circuit characteris-tics, ultra-small nanoscale devices, and for capturing the discrete-dopant-position-induced
fluctuations, a device-circuit coupled simulation approach [5] is employed. The nodal volt-age and loop current in the circuit can be calculated. The formulation of circuit equations is mainly base upon the Kirchhoff’s current law. The circuit nodal equation of 6T SRAM, as illustrated in Fig. 2.8(a), is shown in below:
Node1 : V1 = VDD, (2.13)
Node2 : V2 = VBL, (2.14)
Node3 : V3 = VBL0, (2.15)
Node4 : V4 = 0, (2.16)
Node5 : V5 = VW L, (2.17)
Node6 : Id,P M OS2+ Id,N M OS4= Id,N M OS2, (2.18)
and
Node7 : Id,P M OS1+ Id,N M OS3= Id,N M OS1, (2.19)
The static transfer characteristics of SRAMs are then estimated. Thus, all device and circuit characteristics are obtained without any devices’ equivalent circuit models. The flowchart
2.2 : Simulation Technique 25
for mix-mode simulation method is shown in Fig. 2.9. The characteristics of devices of test circuit are first estimated by solving the device transport equations and using as initial guesses in the device-circuit coupled simulation. The circuit nodal equations of the test circuit are formulated and then directly coupled to the device transport equations (in the form of a large matrix containing the circuit and device equations), which are solved si-multaneously to obtain the circuit characteristics [5]. The flowchart of decoupled PDE is shown in Fig. 2.10. First we solve the nonlinear Poisson equation until it is convergence, and then the current continuity equation of electron and hole is following solved. If the error is less than the tolerance, the program stops and output the initial solution of device’s potential and perform the mixed-mode simulation. We solve the device’s equations cou-pled with circuit nodal equations until it is convergence. Figure 2.11 shows the flow for solving decoupled PDE. First the simulation domain have to be discretized. Applying the finite element approximation to the decoupled PDE, we obtained the nonlinear algebraic equations corresponding to the discretized grid. The Newton linearization method is then used to linearized the nonlinear equations. Finally, the linear algebraic equations can be solved using either direct or iterative method.
Figure 2.4: (a) Discrete dopants randomly distributed in the (80 nm)3 cube with the average concentration of 1.48 × 1018cm−3. There will be 758 dopants within the cube, but dopants may vary from 0 to 14 ( the average number is 6 ) within its 125 sub cubes of 16 nm × 16 nm × 16 nm [(b and (c)].
2.2 : Simulation Technique 27
Figure 2.5: The histogram of the dopants in 125 sub cubes for 16 nm devices. The dopants number can be describe by Gaussian Distribution with a mean of six.
Figure 2.6: The histogram of the dopants in 125 sub cubes for 65 nm devices. The dopants number can be describe by Gaussian Distribution with a mean of 100.
2.2 : Simulation Technique 29
Figure 2.7: The sub-cubes are equivalently mapped into channel region for discrete dopant simulation as shown in MOSFET (a), and the same approach for SOI FinFET (b). The aspect ratio of FinFET is set at 2.0 rather than higher value to examine more critical situation.
Figure 2.8: The schematics of 6T (a) and 8T (b) SRAM cells. Since the bit-line voltage will directly impact the “0” storage node in 6T structure during read operation, additional transistors M1 and M2 in 8T structure is used to avoid the impact thus increase the read stability.
2.2 : Simulation Technique 31
Figure 2.9: The flowchart for the mixed-mode device circuit coupling simulation. The device’s simulation is performed first and get the initial solution of the device’s potential. The
mixed-mode simulation is then executed until the final step.
Figure 2.10: A flowchart of the decoupling algorithm. First we solve the nonlinear Poisson equation until it is convergence, and then the current continuity equation of electron and hole is following solved. If the error is less than the tolerance, the program stops.
2.2 : Simulation Technique 33
Figure 2.11: A flowchart for solving decoupled PDE. First the simulation domain have to be discretized. Applying the finite element method approximation to the decoupled PDE, we obtained the nonlinear algebraic equations corresponding to the discretized grid. The Newton linearization method is then used to linearized the
nonlinear equations. Finally, the linear algebraic equations can be solved using either direct or iterative method.
Intrinsic Parameter Fluctuation on Devices and SRAM Cells
I
n this chapter, the device intrinsic parameter induced characteristics fluctuations on 65-nm to 16-65-nm-gate planar devices were first explored. Then we introduce the stability of an SRAM cell. Finally, the device intrinsic parameter induced SNM fluctuations on 65-nm to 16-nm-gate planar 6T SRAM cells and the correlation between total SNM fluctuation and different transistor pairs were examined.34
3.1 : Physical Characteristics 35
3.1 Physical Characteristics
Figure 3.1(a) shows the ID − VG characteristics fluctuations of the discrete-dopant-fluctuated 16 nm planar MOSFETs, where the solid line shows the nominal case (contin-uously doped channel with 1.48 × 1018 cm−3 doping concentration) and the dashed lines are random-dopant-fluctuated devices. From the random-dopant-number point of view, the equivalent channel doping concentration is increased when the dopant number increases, which substantially alters the threshold voltage as shown in Fig. 3.1(b). The threshold voltage is determined from a current criterion that the drain current larger than 10−7(W/L) ampere. As the number of dopants in channel is increased, the device’s Vth is increased.
The position of random dopants induced different fluctuation of characteristics in spite of the same number of dopants. For the device with the number of discrete dopants, varying from zero to 14, maximum difference of Vth is about 240 mV, which takes an important part in random dopant fluctuation. Furthermore, the magnitude of the spread characteris-tics increases as the number of dopants increases. The physical mechanism can be briefly described by band profile. Figures 3.2(a) and 3.2(b) show the extracted band profile for the nominal case and discrete-dopant-fluctuated case, respectively. For the nominal case, the band profile is smooth. However, there are several potential barriers in the discrete-dopant fluctuated case. These potential barriers are induced by the corresponding discrete-dopants in device’s channel and therefore change the threshold voltage.
Both the randomness of dopants in the channel and source/drain regions may induce the Vth fluctuation in a planar MOSFET device. In order to verify the importance of the RDF in the channel region, the RDF in source and drain regions is also investigated. An example of a 16-nm-gate planar MOSFET with atomistic doping profiles both in the chan-nel and source/drain regions is shown in Fig. 3.3. The actual location of random discrete dopants and the electrostatic potential are illustrated in Fig. 3.3(a). Note that the atom-istic doping in the source/drain of device will not only introduce the electrostatic potential fluctuation but also the variations in the effective length of the channel as shown in Fig.
3.3(b). Fig. 3.4 shows the comparison of RDF effect in channel versus source/drain region.
The Vth fluctuations are 61 mV and 25 mV in channel and source/drain RDF only cases, respectively. Moreover, the fluctuation is only 68 mV for both channel and source/drain RDF cases. It can be seen that the influence of channel RDF on Vth fluctuation is around 90%, which means the Vthfluctuation in a planar MOSFET is dominated by the random-ness of dopants in the channel rather than the source/drain region. Therefore, RDF in the source/drain region can be neglected and won’t be considered in following study.
3.1 : Physical Characteristics 37
Figure 3.1: (a) DC characteristic fluctuations of ID− VGcharacteristics, the solid line shows the nominal case and the dashed lines are random-dopant-fluctuated devices. (b) Vthfluctuated of 16-nm-gate planar MOSFET. The Vthis defined as the gate voltage where the drain current is equal to 0.1 µA. Each symbol indicates one discrete dopant fluctuated case. The threshold voltage may result from discrete-dopant-number and discrete-dopant-position induced fluctuations.
Figure 3.2: Extracted band profile for (a) nominal case and (b)
discrete-dopant-fluctuated case. Several potential barriers in the discrete-dopant fluctuated case are induced by the corresponding dopants in device’s channel and therefore change the threshold voltage.
3.1 : Physical Characteristics 39
Figure 3.3: Example of a 16 nm planar MOSFET with atomistic doping profiles both in the channel and source/drain regions. (a) The actual locations of random discrete dopants and (b) the fluctuation in electrostatic potential are illustrated.
Figure 3.4: Comparison of RDF effect in channel versus source/drain (S/D) region of a planar MOSFET. The Vthfluctuations are 61 mV and 25 mV for channel and S/D RDF only cases, respectively. The Vthfluctuation is 68 mV for both channel and S/D RDF cases. Note that the channel RDF introduced around 90% of Vthfluctuation.