Combined Electromagnetic and Drift Diffusion Models for Microwave Semiconductor Device ()
1. Introduction
Propagation and radiative effects become more and more important for the integrated circuit domain, today submicron semiconductor are operated under high frequencies. This is particularly important for the characterization of interconnected structures loaded at digital and drivers. Indeed, it is well known that electromagnetic compatibility (EMC) and signal integrity (SI) are strongly affected by the geometry of interconnects and by the possibly complex nonlinear/dynamic behavior of the electronic devices collocated at their terminations.
In this paper, the semiconductor analysis is based on the time-domain drift-diffusion model (DDM) in conjunction with electromagnetic model (Maxwell’s equations) [1,2]. Where the active device model (AD model) correspond to the Poisson equation and the carrier transport ones, obtained by splitting the Boltzmann transport equation (BTE) into its first two moments.
In first part of this wok, we develop a modeling and a simulation for a PIN diode.
For this we give a simple discretization scheme for the DDM model using a finite different algorithm to achieve the space-domain integration, whereas Euler backward algorithm is adopted to accomplish time-domain integration. The resulting equations can be written as implicit sparse matrix systems.
In this approach, the three equations are coupled together through the whole process of computation. The well known Gummel’s iteration is most commonly used to solve the coupled problem from the active device model at each time step, whereas each matrix system is solved using an efficient direct UMFPACK method.
In a second part, Finite-difference time-domain (FDTD) method is used to solve electromagnetic equations, where it can be coupled to an efficient solution algorithm by incorporating a “lumped-device” [3]. The design simulation is employed to 3-D structure of the printed circuit including PIN diode (lumped device). The Berenger Perfectly Matched Layer condition (PML) is accomplished for the FDTD computation. The innovation of this condition is that electromagnetic waves of arbitrary incidence, polarization, and frequency are matched at the boundary and subsequently absorbed in the PML layer.
2. Active Device Model
The active device model (AD model) used is based on the moments of Boltzmann’s transport equations obtained by the integration over the momentum space [4]. This model formulates the problem using three dependent variables V, n, and p.
The Poisson’s equation can be employed for the potential within the diode together with the charge conservation for electrons and holes.
(1)
where 
and accounts for the net ionized impurity concentration.
The electron and hole current continuity equation for semiconductor devices can be written as
(2)
(3)
where rSRH is the Shockley-Read-Hall recombination.
(4)
Jn, Jp can be described by a drift-diffusion model:
(5)
(6)
3. Numerical Approach
3.1. Discretization Scheme
The discretization uses a first and second order of (1) in 1-D-finite difference mesh, leads to have
(7)
Euler implicit method seeks to approximate the derivatives in (2) and (3) with regard to the discrete solutions points defined by spatial and temporal cells [5,6]. The electron and hole continuity equation may be discretized in implicit form as follows
(8)
(9)
The variables a1-8 are matrix coefficient which can be found in [Mirzavand], and the indices i describe spatial discretization, k corresponds to the time progression. The mesh size h is limited by the Debye length and ∆t by the dielectric relaxation time The corresponding Matrix systems of (7), (8) and (9) can be written as:
(10)
(11)
(12)
where
A, B and C are the sparse matrix.
3.2. Numerical Whole
In order to derive the iteration procedure, at each time step we relate Gummel’s method to direct UMFPACK method which is well known to converge quadratically [7].
The Equations (10), (11) and (12) are solved in a decoupled manner. Poisson equation is solved at all grid points, followed by electron continuity equation, and then by hole continuity one. Each matrix system is solved using an efficient direct numerical UMFPACK method [8,9].
The whole numerical procedure to calculate the final solution at each time step using Gummel’s iterations can eventually be summarized as follows:
Step 1: Initialize carriers (Evaluate 
and
)
Step 2: Solve (10) (Evaluate
)
For k from 1 to final time iterations Begin Gummel iterations Step 3: Solve (10), (11) and (12) using UMFPACK method Step 4: Test convergence End Gummel iterations Step 5: Set 
and 
End time iterations.
4. Electromagnetic Wave Approximation
The electromagnetic wave propagation can be completely characterized by solving Maxwell’s Equations [10]. These equations are first-order linear the field vectors at any point in the space at any time. The current density distribution is obtained from the solution of the driftdiffusion model [11,12].
(13)
(14)
Here 
account for the contribution of the current flowing along the distributed media, whereas 
includes the contributions of lumped elements. Such a contribution comes from the solution of related device model’s equations.
We consider a locally uniform medium in the vicinity of an x-directed lumped element [13]. A schematic of the coupling between the solvers is given in Figure 1. The entire computational domain is discretized with a FDTD mesh sizes
,
,
along the three Cartesian coordinates [14].
Then the electric field of the lumped element may be defined as [15]:
(15)
where, 
is an impressed current density through which the PIN diode will be incorporated. Such a contribution comes from the time domain solution of the related device model’s equations. The active device is oriented in the x-direction, so the current density is given by:
To preserve the computational efficiency of the “leapfrog” solution scheme, we assume that the device equations are solved at even time-steps; this permit the averaging of current density contribution between the two times steps k and k + 1.
The time step
is bounded by courant stability
criterion [16]
(16)
To simulate the infinite space surrounding the structure, absorbing boundaries conditions are considered. At interconnections, metal surfaces behave, to a good approximation, as perfect electric conductors (PEC). Therefore, we set the tangential component of the electric field to zero on the metal boundaries.
In the cell (i + 1/2, j), we solve drift-diffusion model and Maxwell’s equations, where Equation (15) is adopted to combine between the two model to update the electric and magnetic fields. This algorithm can be efficiently carried out in the tree steps [15,17]:
For k from 0 to final time Step 1: update magnetic field components H using (14)
Step 2: solve semiconductor device using (10), (11), (12) and (15)
Step 3: update electric field components E using (13)
End time iterations
5. Results and Discussions
5.1. Stationary Simulation Results
The PIN diode considered in this simulation is 2.5 µm length. Here, a uniform mesh that covers the 1-D cross section of the diode is used with 500 nodes.
The P and N junctions are doped as Gaussian function with a peak density of 5.1017 At
cm–3 at the two contacts.
The device is biased with a direct potential of Vp = 0.7 v, the DC parameter distributions are obtained by solving the drift-diffusion model.
Figure 2 represents the norm of residual error using Newton’s iterations and direct UMFPACK methods
Figure 2. Norm of residual versus Gummel’s iterations for 1-D p-i-n device, (a) with Newton iteration, (b) with direct UMFPACK method.
respectively. When the UMFPACK method is considered, the CPU time is reduced by 97.34% with a maximum numerical dispersion error of 10 - 11 for this 1D problem. We show clearly that the combined “Gummel-UMFPACK” iteration has significantly better convergence characteristics than the combined “Gummel-Newton raphson” iteration.
Carriers density distributions and potential profile of the device for an applied bias of Vp = 0.7 v is presented in Figure 3.
The static I-V characteristic of PIN diode, as computed for different net doping by means of the numerical model, is shown in Figure 4. Current characteristics are calculated directly by substituting the numerical results in the discretized form of (5) and (6), and then integrate it over the anode contact. The calculated results present a good agreement with those obtained in [15].
5.2. Time Domain Results
At initial time t = 0, the PIN diode is forward biased with the potential Vp = 0.7. As the DC simulation time in-
Figure 3. Sample DC results obtained using the proposed algorithm, (a) potential distribution; (b) Carrier density distribution.
Figure 4. I-V characteristics of the PIN diode, at different doping concentrations.
dependent (the electromagnetic model is not included), where the time step is
.
When a PIN diode is forward biased, holes and electrons are injected from the P and N regions into the intrinsic region. These carriers not recombine immediately. Instead, a finite quantity of charge is stored and results in lowering of the resistivity of the intrinsic region as shown in Figure 5.
Figure 6 represents the variation of the output current versus the time for different values of the potential Vp. We can show clearly that the transient solution converge towards results obtained in the steady state solution (view Figure 4).
5.3. Frequency Domain
The FDTD simulation is performed with uniform grids of space steps ∆x = 0.2 mm, ∆y = 0.2 mm and ∆z = 0.15 mm. The computation domain is divided into a grid of 100 × 50 × 24 and performed for 4000 time steps, It’s located in an unbounded medium, requiring absorbing boundaries placed at a sufficient distance from the circuit to avoid nonphysical reflections.
Figure 7 represents the proposed circuit. The excitation pulse is a Gaussian one with a maximum frequency as 50 GHz. An internal resistor (50 ohm) corresponds at
Figure 5. Electrons distribution versus time.
the 1st port and a matched load (50 ohm) is connected to the 2nd port. Voltage sources and load resistor have been described by analytical models, where the surfaces of the two ports and the ground plane are related by:
Rload = 50 ohm.
The PIN diode is described by numerical solutions obtained in time domain simulation.
Figure 8 shows the special distribution of electric field Ex(x, y, t) beneath the structure at t = 0.882 ps.
The S-parameters of the PIN diode with 50 frequency points is linearly distributed in the frequency—band (1 KHz - 50 GHz). They are calculated and compared with those obtained using finite element software as SILVACOTCAD, which are illustrated in Figures 9 and 10. Time responses for the anode and cathode contacts are calculated using implicit Backward Euler algorithm. A good fit between the proposed model and the finite element simulator (SILVACO-TCAD) is observed.
The insertion loss for the proposed structure is shown in Figure 11, where the real and imaginary parts are presented.
To obtain frequency domain comportment, Fourier Transform algorithm is used in output and input time
Figure 8. Electric field distribution for 250 time steps.
Figure 9. Set of calculated return loss for the PIN diode.
Figure 10. Set of calculated insertion loss for the PIN diode.
Figure 11. Set of calculated insertion loss for proposed structure.
responses, [18] present the details on how such analysis is performed using FFT.
6. Conclusions
In this paper, we have presented an efficient time-domain electromagnetic wave simulator for modeling a PIN diode considered in a micro-wave circuit. The features distributed analysis both the active semiconductor devices and passive interconnection domains.
The mathematical model of the proposed device consists of a set of nonlinear second-order partial differential equations which are solved on uniform mesh using implicit finite difference technique. The descretization of drift-diffusion model is established using backward Euler scheme. The implicit system obtained is solved using the well known Gummel’s scheme related with an efficient Direct UMFPACK method.
Using this approach, the proposed algorithm reaches over 95.34% reduction in CPU time. The modeling results are qualitatively in agreement with theoretical concepts. Ended, the transient solution converges towards the steady state solution which confirms the validity of the proposed algorithm.
In other hand, we have presented an efficient timedomain electromagnetic algorithm for modeling a PIN diode used for mm-wave applications. This model combines both accurate modeling of transport phenomena in PIN diode and electromagnetic wave effects present around the device at high frequencies. The “leapfrog” scheme of 3-D FDTD method is used for the simulation of the propose circuit. Results obtained are in agreement with those obtained by other simulator software like SILVACO-TCAD.
Our conclusion is that, for the following works, we essentially will use the “FDTD” simulator for two essential reasons: the first one: this simulator is also dedicated to other devices analysis. The other reason is we will go on to develop our time domain simulator using CrankNicolson scheme, taking care of the no totally resolved boundary problems, nowadays.
Nomenclature
