{'fraction_non_alphanumeric': 0.051531142560166304, 'fraction_numerical': 0.01870952266730631, 'mean_word_length': 4.2874656520820125, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Cryogenic MOS Transistor Model

This paper presents a physics-based analytical model for the MOS transistor operating continuously from room temperature down to liquid-helium temperature (4.2 K) from depletion to strong inversion and in the linear and saturation regimes. The model is developed relying on the 1D Poisson equation and the drift-diffusion transport mechanism. The validity of the Maxwell-Boltzmann approximation is demonstrated in the limit to zero Kelvin as a result of dopant freeze-out in cryogenic equilibrium. Explicit MOS transistor expressions are then derived including incomplete dopant-ionization, bandgap widening, mobility reduction, and interface charge traps. The temperature dependency of the interface-trapping process explains the discrepancy between the measured value of the subthreshold swing and the thermal limit at deep-cryogenic temperatures. The accuracy of the developed model is validated by experimental results on a commercially available 28-nm bulk CMOS process. The proposed model provides the core expressions for the development of physically-accurate compact models dedicated to low-temperature CMOS circuit simulation.


I. INTRODUCTION
A DVANCED CMOS processes perform increasingly well from room temperature down to deep-cryogenic temperatures (< 10 K) [1]- [4]. At these temperatures the ideal switch with a step-like subthreshold slope comes within reach [5]. Furthermore, cryo-electronics [6]- [8] can provide an interface with superconducting devices on the quest for exascale supercomputing [9]. Ultimately, quantum-engineered devices controlled by cryo-CMOS circuits can bring new functionality to existing computing technologies [10], [11].
Large-scale integration of silicon spin qubits [12], [13] and cryo-CMOS control circuits is envisioned to take solid-state quantum computing to the next level [14]. Digital, analog, and RF CMOS circuits [15]- [17] are then required to operate at millikelvin temperatures for initialization, manipulation, and read-out of the qubits, as well as error correction [18], [19]. Since the cooling power at millikelvin temperatures is reduced, the system could feature a cryogenic temperature gradient, where the control circuits operate at a higher cryogenic temperature than the qubits, e.g. 4.2 K [15]. However, the optimal design of power-hungry and thermal-noise dissipating circuits operating in close proximity to the qubits is yet to be explored. In this context, the main hurdle to overcome is the lack of compact MOS transistor models in circuit simulators, remaining physically accurate below 10 K [15], [17].

II. CRYO-MOS TRANSISTOR MODELING
The low-temperature circuits developed for spacecraft [20], [21], scientific equipment [22], ultra-low-noise detectors [23], cryobiology [24], and others, have been custom-designed relying on a semi-empirical approach. This approach requires laborious and expensive low-temperature measurements to extract model parameters for tuning room-temperature compact models to the target low-temperature [23], [25], [26]. Empirical temperature-scaling laws have been added to the room temperature physics-based MOS transistor model [27], [28] to capture cryogenic operation down to 4.2 K [29]- [31]. However, the discrepancy between the measured value of the subthreshold swing for a long device at 4.2 K (≈ 10 mV/decade) [3], [4], [32], and the theoretical thermal limit, U T ln 10 (≈ 0.8 mV/decade) reveals that something more fundamental is missing. As we will demonstrate along this paper, important physical phenomena at low temperatures such as interface trapping [28], [33] and incomplete ionization [34], [35] have not been properly included to date. Furthermore, the intrinsic carrier concentration, n i , takes on extremely small values below 10 K, causing arithmetic underflow in implemented analytical expressions or convergence problems in computer-aided-design simulations [36]- [38]. Therefore, standard references on semiconductor devices treat only the cryogenic equilibrium condition in bulk semiconductors above 10 K [27], [28], [39]. Analytical device-physics models, starting from the Poisson equation at low temperature, leave a gap unfilled between the zero-Kelvin approximation and 77 K [40]- [43].
In this work, we develop a MOS transistor model valid from room temperature (RT) down to deep-cryogenic temperatures, entirely based on physics principles and validated with experimental results. We start by verifying the continued validity of the Boltzmann statistics down to the deep-cryogenic regime.

III. MOS ELECTROSTATICS FROM RT TO 4.2 K
We model a long, planar n-channel MOS field-effect transistor in silicon, depicted in Fig. 1. Uniform operation across the width of the transistor is assumed and the gradual channel approximation is adopted. The electrostatics can then be described by the 1D Poisson equation [27], [28].
where q is the elementary charge, ε si the silicon permittivity, and ψ (E F − E i )/q the potential, with E F the Fermi-level and E i the intrinsic energy level. The first term on the RHS of equation (1) represents the electron contribution, n, the second term the hole contribution, p, and the third term the ionized dopant-contribution, N − A . 1) Incompletely-ionized dopants: Under thermal equilibrium, both at room and cryogenic temperatures, the majority carrier concentration can defer from the implanted doping value, N A , due to incomplete ionization of the dopants. In cryogenic equilibrium, incomplete ionization is strong and known as freeze-out, since thermal dopant-ionization is very low [35]. However, during MOS operation, also field-assisted ionization comes into play. Fermi-Dirac statistics provides a fundamental way to model incomplete ionization which includes both dopant-ionization mechanisms. The concentration of ionized dopants, N − A , is then equal to the total concentration of implanted dopants times the Fermi-Dirac occupation probability of the acceptor energy where the electron quasi-Fermi-level is given by E F,n = E F −qV ch . The RHS of (2) is obtained by replacing E A −E F,n with E A − E i + E i − E F,n in the exponential term, and by defining an acceptor potential, ψ A (E A − E i )/q, as indicated in Fig. 1. The channel voltage, V ch , denotes the shift of the quasi-Fermi-potential due to the drain-to-source voltage, V DS . The second expression in (2) highlights the two dopant-ionization contributions, i.e. the potential (fieldassisted ionization [35]) and temperature (thermal ionization). The acceptor-site degeneracy factor, g A , is set to four due to fourfold degeneracy (heavy-light hole, spin up-down) [28], [39]. Note that setting g A to zero is equivalent to assuming complete ionization.
2) Mobile carrier concentrations: Since n and p given by Fermi-Dirac statistics in (1) require numerical integration over energy, this inhibits explicit solutions for the charge densities and current in the MOS transistor. Expressing n and p using Boltzmann statistics allows to obtain such relations. However, the validity of the Maxwell-Boltzmann approximation down to deep-cryogenic temperatures is questionable. It has been reported [38], [42] that semiconductors become strongly degenerate at deep-cryogenic temperatures, preventing its use. This is however inconsistent with the zero-Kelvin limits of the Fermi-level position in the bandgap derived by Pierret [39]. Therefore, in the next subsection we aim to verify the Maxwell-Boltzmann approximation down to deep-cryogenic temperatures.
3) Verification of Boltzmann statistics: We numerically calculate the position of the equilibrium Fermi-level, E F , down to 100 mK relying on Fermi-Dirac statistics in an extrinsic bulk  Fig. 2a), leading to bulk freeze-out according to (2). When E A bends under E F near the surface, the acceptor dopants become rapidly completely ionized due to field-assisted ionization. The quasi-Fermi potential is not taken into account in this figure. semiconductor, e.g. p-type silicon. In this case, the Poisson equation imposes the charge neutrality, p p = N − A , where p p is expressed by Fermi-Dirac statistics [28], [39] and N − A by (2). This yields an implicit equation for E F , which is solved numerically at each temperature and doping value using an extension of the arithmetic precision. As illustrated in Fig. 2a, below 120 K, E F remains off the valence band edge with an offset larger than 3kT for doping values below the degenerate limit (i.e. N A = 4 × 10 18 cm −3 in Si:B) [28], [39], [45]. Note that this is predicted correctly only when incomplete ionization is taken into account. Complete ionization (g A = 0) would predict an offset smaller than 3kT for N A = 10 18 cm −3 , and hence a degenerate semiconductor. It should therefore be emphasized that incomplete ionization maintains the non-degeneracy of a highly-doped semiconductor at temperatures down to 100 mK. Furthermore, near zero Kelvin, E F tends to saturate at (E A − E v )/2 for all considered doping values. This corresponds to the zero-Kelvin limit in Pierret [39] assuming Boltzmann statistics. Using the now validated Maxwell-Boltzmann description for Considering the temperature dependency of N v [28], [39], while taking the limit of (3) to 0 K, leads to lim Performing the same numerical E F -calculation for an intrinsic semiconductor, the extremely small value of n i can be verified relying on Fermi-Dirac statistics. The Poisson equa- tion then imposes the charge neutrality, n = p = n i , where n and p are given by Fermi-Dirac statistics. As illustrated in Fig. 3, this yields n i -values lying outside the range of IEEE double-precision arithmetic (10 −308 − 10 308 ), e.g. at 4.2 K, ≈ 10 −678 cm −3 . Therefore, an extension of the arithmetic precision will also be used in the remainder of this work based on Boltzmann statistics, since the carrier concentrations are then expressed through n i .

B. Poisson-Boltzmann equation
Using the Maxwell-Boltzmann approximation of n and p, validated down to deep-cryogenic temperatures in the previous section, we combine the 1D Poisson equation with Boltzmann statistics, which leads to where U T kT /q is the thermal voltage. The first term on the RHS of equation (4) represents the electron contribution, n, and the second term the hole contribution, p. The intrinsic carrier concentration is given by , with E g the bandgap, and N c and N v the effective density-ofstates in the conduction band and valence band respectively. The temperature-dependency of E g as described by Varshni [44] is used. The extremely small, but finite value of n i at deep-cryogenic temperatures cannot be assumed zero-which would be equivalent to the zero-Kelvin approximation [40] or considering f (E) as a step function-since this leads to zero mobile carrier concentrations independently of the potential. This is irreconcilable with the observed field-effect and correct functioning of the MOS transistor at 4.2 K [3], [5].
For U T small, the exponential factor has a very big dynamic range when ψ changes during MOS transistor operation, large enough to overrule n i in the multiplication.
Integrating (5) from bulk to surface with E = −∂ψ/∂y and In (6) the additional potential dependence due to field-assisted ionization of the dopants can be straightforwardly integrated as well, i.e. by replacing in the numerator of the third term and splitting the resulting integral. This gives an expression for the square of the electric field at the surface, where ψ b (E F,b − E i )/q is the bulk potential and ψ s (E F,s −E i )/q the surface potential, as indicated in Fig.1. E F,s denotes the Fermi-level at the surface, and E F,b the Fermilevel in the bulk. The logarithmic term in (7) is the contribution the Fermi-Dirac ionization probability at the surface, and the Fermi-Dirac ionization probability in the bulk, assuming that V ch is zero in the bulk. Both ionization probabilities are qualitatively shown in Fig. 1. If complete ionization is assumed, then f s (E A ) = f b (E A ) = 1 and the incomplete ionization term cancels in (7), leading to the expression widely-used at RT [27], [28]. The surface-ionization probability f s (E A ) is plotted in Fig. 4 as a function of thermal and field-assisted ionization. Immediately evident is that freezeout at the surface (arbitrarily defined when f s (E A ) < 0.2) is only present when the temperature is below ≈ 50 K and the potential is close to the flatband condition (ψ s ≈ ψ b ). Above ψ b , the ionization probability rapidly transitions to one due to field-assisted ionization. This transition corresponds to the bending of E A under E F at the surface in Fig. 1. Therefore, complete ionization is a valid approximation even at deep-cryogenic temperatures, although the shift in E F due to incomplete ionization (Fig. 2a) should be taken into account since it affects the threshold voltage. This E F -shift can be quantified by using f b (E A ) from (9) in the bulk charge neutrality condition, p p = N − A , which leads to the quadratic equation The second term in (10) is the shift of E F by including incomplete ionization, which is only dependent on temperature and doping. Assuming complete ionization, i.e. g A = 0, the well-known U T ln(n i /N A ) is obtained.
2) Derivation of the charge densities: Applying the Gaussian law over the semiconductor body in Fig. 1, the total semiconductor charge density per unit area, Q sc , is obtained by Q sc = −ε si E s , with E s given by (7). The obtained Q sc is plotted in Fig. 5a at RT, 77 K, and 4.2 K. For 77 K and 4.2 K, small kinks are noticeable close to ψ b due to the transition from incomplete to complete ionization when E A bends under E F at the surface (E F,s ), or equivalently, ψ s becomes less negative than ψ A . Above this transition, f s (E A ) ≈ 1 according to (2). At RT, E F lies above E A in the flatband condition (see Fig.2b) and hence no transitional kink is noticeable. There is however a ψ b -shift also at RT due to incomplete ionization according to (10). Note that for complete ionization (dashed lines) no kinks are observed since the logarithmic term cancels in (7). Assuming the charge-sheet and fully-depletion approximations [27], the fixed charge density per unit area, Q f , is given by Relying on the charge neutrality, the mobile charge density per unit area, Q m , can be obtained from Q m = Q sc − Q f , resulting in Eq. (12). Q m is plotted in Fig. 6a for RT, 77 K, and 4.2 K. As can be observed in this figure, incomplete ionization does not affect the turn-on rate of Q m , but contributes a small decrease in the charge threshold voltage, due to the shift of the E F -position closer towards the conduction-band edge as shown in Fig. 2a and derived in (10). Therefore, from this section we conclude that incomplete ionization cannot explain the offset between the measured subthreshold swing at 4.2 K and the thermal limit. As we will show in the next section, the temperature-dependent occupation of interface charge traps can degrade the subthreshold swing down to 4.2 K.
3) Interface charge traps: Defects and lattice breaking at the oxide-semiconductor interface introduce trap energy levels, E t , in the bandgap which degrade the control of the gateto-bulk voltage, V GB , over the channel. In what follows, the Fermi-Dirac occupation of interface traps, f (E t ), is included in the surface-boundary condition and the effect on the Q m turn-on rate is analyzed at 4.2 K. The surface-boundary condition, i.e. the link between V GB and ψ s , is given by where C ox is the oxide capacitance per unit area, and V FB is the flatband voltage, given by V FB φ ms − Q it /C ox [27], [28]. Here Q it is the interface trap charge density per unit area. We consider the summation of the discrete acceptor trap energy levels (all donor states are occupied and neutral during turn-on in nMOS [28]). Each discrete trap-energy-level, E t,j , at position j in the bandgap has its particular N it,j -value assigned to it, where N it is the density-of-interface-traps per unit area. Q it can then be expressed as Q it = −q N j N it,j f s (E t,j ) where N is the number of interface traps, and is the Fermi-Dirac occupation probability of the trap energy level E t,j . The RHS of Eq. (13) is obtained by defining the trap potentials, ψ t,j (E t,j − E i )/q [46]- [48]. This leads to the flatband voltage Plotting Q m from (12) versus V GB at 4.2 K in Fig. 6b, including four interface traps close to the conduction band, reveals how each interface trap degrades the turn-on of Q m separately, as well as the combined effect of the sum of the interface traps.

IV. CURRENT DERIVATION
To derive the current in the linear regime, this core-model assumes drift-diffusion transport, and does not include ballistic nor quantum transport. To verify the drift-diffusion transportmechanism at cryogenic temperatures, the proposed model for the drain-to-source current will be experimentally validated in Section V. Neglecting the hole-contribution to the current, the expression for the total drain-source current is given by where the electron mobility µ n is assumed constant along the channel, and W/L is the device aspect-ratio, as indicated in Fig. 1. In the linear regime, Q m can be assumed independent of V ch . In this case, the total drain-source current is given by I DS = −µ n (W/L)Q m V DS . In saturation, the integral over V ch cannot be readily solved. Fig. 4. Dopant ionization at the surface is an interplay between thermal ionization (T ) and field-assisted ionization (ψs − ψ b ). Freeze-out is assumed when 20 % of the dopants are ionized. This happens only when T is below ≈ 50 K and close to the flatband condition (ψs ≈ ψ b ). When ψs increases, a rapid transition takes place to complete ionization for all temperatures. In the flatband condition (ψs = ψ b ) the ionization probability at the surface is only due to thermal ionization and equals the ionization probability in the bulk, f b (E A ).

Therefore, starting from the drift-diffusion equation gives
Assuming a linearization of the mobile charge density with respect to the surface potential at constant gate voltage [49] in (15), i.e., Q m = mC ox (ψ s −ψ P ), with m ∂(Q m /C ox )/∂ψ s and ψ P the pinch-off potential, and integrating, results in an expression for the total drain-source current in saturation, Q m,S and Q m,D are obtained from (12), setting V ch to zero and V DS respectively. At cryogenic temperature an improvement in the low-field mobility, µ 0 , is observed due to a reduction of the phonon scattering [1], [6]. In addition, the mobility reduces at higher gate voltage due to surface roughness scattering at high vertical electric field [1]. This mobility reduction can be modeled by µ n = µ 0 /(1 + θV GB ) where θ is the mobility reduction factor.

V. EXPERIMENTAL RESULTS AND DISCUSSION
Room temperature and cryogenic measurements were performed on devices fabricated in a 28-nm bulk CMOS process. The full set of measurements, measurement set-up, and characterization were previously reported in [3], [4]. After measuring at RT, the samples were immersed into liquid helium (4.2 K) and liquid nitrogen (77 K) baths with a dipstick. Fig. 7a favorably compares the model with the linear transfer characteristics (V DB = 20 mV) measured at RT and 4.2 K on a long nMOS device with W/L = 3µm / 1µm, in linear and logarithmic scales. The extracted µ 0 -values from the model are in accordance with the characterization performed in [3]. Furthermore Fig. 7b analyzes the effect of incomplete ionization, interface traps, and mobility, on the current at 4.2 K. Note that incomplete ionization reduces the threshold voltage, and interface traps can degrade the subthreshold swing (SS) to ≈ 10 mV/decade. The strong increase in the mobility increases the on-state current at 4.2 K. Fig. 8 validates the model for the current in saturation (|V DB | = 0.9 V) using the measurements performed at RT, 77, and 4.2 K on a long pMOS device with W/L = 3µm / 1µm, in linear and logarithmic scales. The metal-semiconductor work function difference, φ ms , increases in absolute value at lower temperatures according to the change in E F -position (Fig. 2).

VI. SUBTHRESHOLD-SWING DERIVATION
In this section a SS-expression including incomplete ionization and temperature-dependent interface trapping is derived. Incomplete ionization is included to prove the minimal influence on SS shown in Fig.7b. The temperature-dependency  , green), and liquid-helium temperature (4.2 K, blue) including incomplete ionization (solid lines) or assuming complete ionization (dashed lines). The potential is swept starting from the bulk potential, ψ b , calculated at a given temperature and doping according to (10). Horizontal arrows show the shifts in ψ b by including incomplete ionization at a given temperature. Skewed arrows in the insets indicate the kinks at 77 and 4.2 K due to the transition from incomplete to complete ionization when E A bends under E F (Fig. 1). of interface-trap occupation, f s (E t ), allows to obtain the ∆ SS -offset of ≈ 10 mV/dec above the thermal limit, U T ln 10, previously observed on long-channel devices [3], [4]. The subthreshold swing, SS, is usually expressed as nU T ln 10, where the non-ideality factor or slope factor, n, is given by (∂V GB /∂ψ s ), describing the deviation from the thermal limit. Assuming f s (E t ) from (13) to be one, (∂V GB /∂ψ s ) [27], [50]. However, at 4.2 K, and assuming the highest possible doping value below the degenerate limit, a large N it -value in the order of 10 13 cm −2 is extracted to accommodate for a SS of ≈ 10 mV/decade [33], [51], since N it becomes multiplied with U T in this expression. However, it should be emphasized that in the used expression for SS, the temperature dependency of interface-trap occupation is not taken into account. Relying on drift-diffusion transport in the linear regime and assuming µ independent of V GB , the subthreshold slope, SS −1 , is given by The factor (1/Q m )(∂Q m /∂ψ s ) is found from (12) by considering Q m Q f or Q sc ≈ Q f in the subthreshold region. After some mathematical manipulation, we find Merging (11) and (18), and inverting, gives The following relation can be derived for (a) in the above equation (see Appendix),  Plugging this in (19) one finds where follows from the surface-boundary condition derived in Sec. III-B.3. In the subthreshold region, far above the flatband condition, f s (E A ) = 1 can be assumed (Fig. 4), and U T 2(ψ s − ψ b ), leading to SS = U T ln 10(∂V GB /∂ψ s ) with Taking the derivative of (13), (23) becomes Note the appearance of a factor 1/U T in the third term on the RHS of (24). Placing a discrete interface trap, ψ t,j , at each ψ s -value in the subthreshold region, and assuming a uniform N it,j -value for each trap, leads to The first two terms in (25) yield the non-ideality or slope factor without interface traps, n 0 . The SS-expression becomes SS = n 0 U T ln 10 + qN it C ox g t (1 + g t ) 2 ln 10 (26) where the second term on the RHS is the sought ∆ SS -offset. The non-ideality factor n 0 has an upper bound of 2 mainly related to doping. Therefore, at 4.2 K, the first term on the RHS of (26) is limited by 1.6 mV/decade. Note that N it does not become multiplied with U T in (26). Therefore, assuming a reasonable value for N it = 3 × 10 11 cm −2 , the second term gives ≈ 9 mV/decade (with C ox = 20 mF m −2 and g t =4).
Together they yield the SS-degradation observed on a long nMOS device at 4.2 K. At 77 K, a similar calculation using n 0 = 1.08 and N it = 3 × 10 11 cm −2 gives 25 mV/dec, corresponding to the subthreshold swing measured [3] on a long pMOS device at 77 K (Fig. 8).

VII. CONCLUSION
A theoretical MOS transistor model is developed valid from room temperature down to liquid-helium temperature. The model relies on Boltzmann statistics, verified in the limit to zero Kelvin, and includes incomplete ionization, interface traps, bandgap temperature dependency, and mobility reduction. It is evidenced that incomplete ionization maintains the non-degeneracy of a semiconductor at deep-cryogenic temperatures, and leads to a decrease in the threshold voltage on top of the overall increase due to Fermi-Dirac distribution scaling. The Fermi-Dirac temperature-dependency of interface-trap occupation degrades the subthreshold swing down to 4.2 K. An expression for the subthreshold swing including incomplete ionization and temperature-dependent interface trapping is derived. The proposed model builds the indispensable physical foundation for future low-temperature CMOS circuit design.

APPENDIX
Starting from Q m +Q f = Q sc , we can write (Q m +Q f ) 2 = ε 2 si E 2 s , with E s given by (7). Solving a quadratic equation for Q m leads to where we neglected the exponential term in ψ b /U T . In the subthreshold region (Q m Q f ), √ 1 + x can be approximated by 1 + x/2 for x → 0. Using (11) for Q 2 f leads to (20).