## 1. Introduction

Natural radiation at ground level, including both terrestrial cosmic rays and telluric radioactivity (alpha decay from radioactive ultra-traces in materials), is considered today as a major reliability issue for integrated circuits (IC), since they are increasingly sensitive to this radiation as long as CMOS technologies scale down [1]. The basic mechanism taking place when an ionizing particle crosses a circuit is the generation in the semiconductor region of a high-density charge track; the generated charge can be collected at circuit level through biased contacts and especially reverse-biased drain junctions and create parasitic transient currents at the circuit nodes. In the case of a memory cell, the collected charge may be sufficient to induce the cell upset; this phenomenon is called single event upset (SEU). When combinational logic is concerned, the charge induced by the particle may lead to a single event transient (SET) [2].

Researchers have always used modeling and simulation approaches to better understand the physical processes and to predict the consequences on the circuit operation of such single events. The study of these phenomena, and more precisely, the computation of the electrical response of devices and circuits submitted to ionizing radiation imply the accurate modeling of various physical mechanisms at different scales, as charge generation, transport and collection of charges within the circuit. The solutions developed to solve this kind of problem include full numerical methods (such as TCAD) and approximate analytical models used in circuit simulators [3]. Full numerical methods, such as TCAD, exactly solve this coupled problem (radiation + electrostatics + transport) and offer a very accurate solution. But the huge computation time and computer resources needed to solve this high complexity problem restrict its application to small simulation domains, limited to, at most, several devices (very small circuits). In order to overstep these limitations, several approximated solutions have been developed in the literature: these approaches generally calculate the current transients induced by single events (SETs) in devices and circuits by decoupling the radiation-induced charge transport from the electrostatic problem. The so-called diffusion-collection model is one of the most interesting and physically accurate approaches developed for this purpose [4–6]. This model makes the assumption that the main carrier driving mechanisms is a pure spherical ambipolar diffusion in the semiconductor region. Moreover, the model can be easily implemented since an ionizing particle track can be divided into a series of discrete punctual charges to compute (or to analytically derive) single or multiple node charge collection. Although interesting improvements have been made [6, 7] and in spite of its capability to be parallelized [8], the diffusion collection model does not correctly take into account the electric field effects in the vicinity of the collecting contact, which is a serious intrinsic limitation. To overcome this drawback, the diffusion coefficient D (i.e. the mobility µ) and the carrier collection velocity v are generally fine-tuned on results obtained from numerical simulation (TCAD) in order to fit the transients for a series of typical radiation events. But, with this procedure the value of D (or µ) does not quantify the real diffusion process (the diffusion will be overestimated or underestimated), since, it artificially reproduces the combination of both the diffusion process and the collection by a reverse-biased junction. Another limitation of this approach is its inability to take into account a non-constant electric field or the nodes bias changes during the collection process.

This chapter presents a new computational method for accurately modeling single-event transient current and charge collection at circuit level. This approach is called random-walk drift-diffusion (RWDD) [9, 10] and consists in a fast Monte Carlo particle method based on a random-walk process [11] that includes both the diffusion and the drift of carriers in an electric field that is nonconstant in both space and time. This method has been successfully used in previous works, as shown in Refs. [12, 13]. In the present work, we considerably improve this approach by a series of developments that make it capable of fast and intensive simulation of full circuits with an accuracy degree similar to that of the TCAD simulation. The considered improvements are briefly described below, but they will be extensively detailed in the following sections. Firstly, the method models the charge carrier transport and collection process using exclusively physical parameters (for example for the carrier mobility or the minority carrier lifetime) and physical constants. In particular, transport and collection processes are modeled by the same mathematical equations, without any fitting parameter. Secondly, the new method implements a model fully derived in C++ environment, which authorizes the use of powerful capabilities such as advanced structures (classes and containers) to model the different problems to solve (geometry, particle track, charge transport, collection, recombination and extraction). Finally, the third important improvement concerns the intrinsic ability of the model to be parallelized on graphic unit processors (GPUs) for charge carrier transport and dynamically coupled with SPICE to take into account the impact of charge collection on biased nodes not considered as stand-alone contacts but embedded in a circuit and connected to other nodes.

This chapter describes in-depth the RWDD approach and illustrates all of these points in the three following sections. Section 2 details the physical insights of the model. The practical implementation of the method using an object-oriented programming language and its parallelization on GPUs is presented in Section 3. Besides, the capability of the approach to treat multiple-node charge collection is presented. Section 4 details the coupling of the model either with an internal routine or with SPICE program for circuit solving. Finally, Section 5 illustrates this proposed approach at device and circuit-level, considering four different test vehicles in 65 nm technologies: a stand-alone transistor, a CMOS inverter, a SRAM cell and a flip-flop circuit. RWDD results are compared with data obtained from a full three-dimensional (3D) numerical approach (TCAD simulations) at transistor level. The importance of the circuit feedback on the charge-collection process is also demonstrated for devices connected to others.

## 2. Random-walk drift-diffusion (RWDD) model

This section introduces the basis of the RWDD model and the modeling methodology used to describe the main simulation steps, i.e. the charge deposition induced by the passage of an ionizing particle at silicon level, the radiation-induced charge transport within the structure and the computation of the collection current on the collecting node(s).

### 2.1. Charge deposition

In the RWDD approach, an ionizing particle track crossing a circuit at silicon-level is modeled as a series of charge packets (**Figure 1**) spread along a straight segment whose length equals the ionizing particle range (R) in target material. The linear density of the charge packet along the particle track takes into account the non-constant linear energy transfer (LET) of the particle. SRIM tables are used to compute both LET and range using appropriate numerical functions [14]. The accuracy of this charge discretization is ensured by the degree of “granularity” [11] that can be fine-tuned by selecting the packets size, in practice from 1 to 100 elementary charges.

### 2.2. Charge transport modeling

The transport of each charge packet starts immediately after the particle crosses the device (**Figure 1**). For modeling the charge transport, the RWDD approach implements the popular drift-diffusion model [15], where the electron (**J**_{n}) and hole (**J**_{p}) current densities are computed as the sum of two components: the drift current component (first term) that is driven by the electric field and the diffusion component (second term) that models the current induced by the gradient of carrier concentration. **J**_{n} and **J**_{p} are given by:

In these equations, n and p are, respectively, the electron and hole densities; µ_{n} and µ_{p} are, respectively, the electron and hole mobilities; E is the electric field; and D_{n} and D_{p} are, respectively, the diffusion coefficients that may be calculated from the carrier mobility using the Einstein’s equation:

In Eq. (3), k_{B} is the Boltzmann constant and T is the carrier temperature; as the carrier gas in the drift-diffusion approximation is assumed to be in thermal equilibrium, T is equal to the lattice temperature.

The next calculation step in conventional full numerical methods (TCAD) is to inject the current densities given by Eqs. (1) and (2) into the conservation laws for electrons and holes also called continuity equations; they are next self-consistently solved with Poisson’s equation. For this purpose, the simulation domain is meshed, and all previous equations are discretized on this mesh grid and then solved. Contrary to this procedure, in the RWDD model, no meshing is need since charge packets have continuous coordinates. Also in RWDD, a random-walk algorithm [11] is used to model the diffusion process, and the drift-induced current is directly calculated using the electrical field present in the considered region. For a charge packet situated at the position (x, y, z) at time t, its new position at time t+δt is given by (x+δx, y+δy and z+δz), where dx, dy and dz are calculated as follows:

In Eq. (4), N_{1}, N_{2} and N_{3} are three independent normal random numbers, D is the diffusion coefficient, µ is the carrier mobility and E(E_{x}, E_{y}, E_{z}) is the electric field vector at the corresponding position and time.

### 2.3. Collection current computation

At each time step of the simulation, the radiation-induced collection current is computed from the transport dynamic of minority charge carriers described in the previous Section 2.2. For the estimation of this current, two main procedures may be employed: the first technique is to use the semi-conductor transport equations and the second one is to employ the Schockely-Ramo’s theorem. Simulation tools used in microelectronics generally consider the first option; the transport equations are implemented in TCAD simulator that numerically solve them self-consistently with Poisson’s equation. This approach considers the free-charge carrier distributions as continuous functions in time and space coordinates. The second option is generally used in instrumentation or high-energy physics for the calculation of detector responses to radiation events. In this second approach, the Ramo’s theorem is used to treat each carrier considered individually and all the interesting effects due to particular carriers are summed [16]. In our work, we implemented the first formalism (transport equations) by applying the continuity equation at the collecting (drain) contact. The transient current at the collecting node is directly computed from the number of carriers Δn that reach this contact during the time step Δt, i.e.

In this expression, the displacement current is neglected; that is a reasonable approximation in this case [17]. This collection current is then injected in the electrical simulation to model the circuit response.

### 2.4. Model limitations

Although the model embeds a microscopic physical description of charge generation, transport and collection, the present approach suffers from two evident limitations that should be amended in a future enhanced code version: (1) the calculation is not self-consistent with the electrostatic potential (i.e. the electric field) and (2) it does not implement electrostatic interactions taking place between charge packets. For this reason, the RWDD model should especially fail in high carrier injection conditions, i.e. at high LETs. These drawbacks will be carefully estimated and quantified in a future work. Another interesting point should be to compare the two methods for computing the collection current to also evaluate the impact of the displacement current on transient characteristics.

## 3. Model implementation

### 3.1. C++ implementation

To implement the RWDD model described previously, we choose a C++ programming environment, an object-oriented programming language that offers considerable advantages in terms of advanced structures, such as classes, objects and containers. The C++ code implementing the RWDD approach was named RWDDCPP.

In this code, the particle track is defined as a C++ container including all the charge packets described as independent objects. The members of the charge packet class include the geometrical coordinates of the packet and the electrical charge amount per packet. The container content may change during the simulation due to two different mechanisms: (1) minority carrier recombination or (2) carrier extraction that correspond to particles that escape the simulation domain. For modeling the carrier recombination mechanisms, we use a simple exponential law that adjust the number of charges as a function of time: N(t) = N_{0} exp(−t/τ). In this exponential law, N_{0} is the initial number of charges deposited by the particle at t = 0 s and the time constant τ is equal to the carrier lifetime.

In RWDDCPP code, the circuit geometry is implemented as a series of 3D rectangular boxes that represent respectively the substrate, the source and drain contacts at the silicon level and the different wells; all these elements constitute the Front-End-Of-Line (FEOL) structure. A simplified Back-End-Of-Line (BEOL) structure may be also modeled, as a stack of insulating material and metal layers. The fine modeling of the reversely biased drain contacts collecting the minority carriers created along the particle track is the most important improvement of the RWDD model at the circuit level. A special “drain class” has been developed, embedding in the same C++ object both the drain contact geometry and the 3D distribution of the electric field induced by the drain bias. At this level, different doping profiles can be taken into account for the p-n junctions: abrupt, gradual or user-defined profiles.

An object of the class “drain” is shown in **Figure 2**; this class is defined by its doping profile (an abrupt junction in this case), its geometrical dimensions (W, L, h), the doping levels of the p and n+ regions (N_{A} and N_{D}, respectively) and the bias (V_{R}) applied to the collecting contact. A 3D analytical model is used to compute the space charge region, since an abrupt p-n junction is considered in this example. The electric field (E) developed in this space charge region and its width (W_{SCR}) are functions of the drain bias potential and can be calculated using the following relations:

where V_{bi} is the internal potential of the junction and ε_{Si }is the permittivity of silicon.

It must be noted that both W_{SCR} and E can dynamically vary during the transient simulation as a function of the effect of carrier collection on the circuit node potential V_{R}, as illustrated in the following.

### 3.2. GPU parallelization

As explained in the previous section, in the RWDD approach, the behavior of each packet of charge is computed independently of the other charges. These processes being independent, the calculation task of the charge transport can be easily parallelized on a graphic processing unit (GPU) whose internal architecture is perfectly adapted to such a massive parallelism. Moreover, RWDD model needs to implement a random number generation procedure that is usually very time-consuming; if the random numbers are independent, this task can also be easily parallelized on GPU. These two characteristics of the RWDD model are expected to offer a considerable computation speed, since the number of parallelizable tasks in RWDD is relatively high. In this work, we used the CUDA programming framework proposed by NVDIA (CUDA version 5.5 64 b) [18]. Random numbers are generated using the CUDA cuRAND random number generation library (Mersenne Twister) within the parallel kernels running on the GPU. **Figure 3** illustrates this parallel implementation of the RWDDCPP code. For the purpose of this work, we separately tested the RWDD algorithm implemented in series on a CPU and in parallel on a GPU. Tests have been performed on a machine equipped with a 3.2 GHz Intel Core i5 and with a NVDIA GeForce GTX 675 MX (960 threads at 600 MHz). Our results show that for the simulation of 100,000 charge packets with 1000 time steps the speedup on GPU was 140× with respect to the classical CPU implementation.

### 3.3. Multiple node charge collection

Since the RWDD model has been implemented in an object-oriented language (C++) using dynamic containers, as previously explained in 3.1, complex circuits with an arbitrary number of sensitive areas and collecting nodes can be simulated, intrinsically considering, in this case, multiple node charge collection in the simulation process. The simulation begins with the initialization of the structure and the definition of the circuit nodes corresponding to the different simulated drains. In a second step, the electrical potential of each node is extracted from a steady-state circuit simulation. The values of these potentials are used to initialize the electrical field in the complete structure. Starting the time-domain analysis, at each time step of the transient simulation, the magnitudes of the current sources corresponding to the different collection processes (simultaneously occurring at different circuit nodes) are updated as a function of the number of collected charge packets. **Figure 4** illustrates this multiple node charge collection process in the simple case of two adjacent drain junctions. The time evolution of the charge packets induced by a horizontal ionizing particle in this two-node structure is shown at five different times after the ionizing particle impact. As evidenced in **Figure 4**, the collection efficiency of the two drains differs because they correspond to two different nodes in the circuit and their electrical potential (i.e., internal electric field) is not the same in this example.

## 4. Circuit solving

Once the current transients have been evaluated using the RWDD model on the different collection nodes of interest in a given device or circuit, one must evaluate the transient electrical response of the circuit. This latter can be directly computed using an internal routine (solving the Kirchoff’s circuit laws) or with an external circuit simulator program. We examine here these two solutions, both implemented in our code RWDDCPP.

### 4.1. Internal routine

For circuit architectures composed of only a few connected devices, such as CMOS inverters or SRAM cells, the steady-state circuit electrical solution can be simply solved considering Kirchoff’s circuit laws and compact models for transistors. In the following, the method is illustrated for the solving of an inverter and a SRAM cell. **Figure 5** (left) defines the different terminal voltages and the source-to-drain current for the NMOS and PMOS transistors; **Figure 5** (right) shows the circuit schematic of a SRAM cell consisting of two cross-coupled CMOS inverters. Considering that the particle strikes the OFF-state NMOS transistor termed NMOS_{1} (initial conditions V_{1} = 0, V_{2} = V_{DD}) the time variations of potential V_{2} for the sole and isolated inverter #1 can be written as:

where I(t) is the current pulse given by Eq. (5) due to the particle strike and collected on the node 2, I_{DN1} and I_{DP1} are the currents of the NMOS_{1} and PMOS_{1} transistors, respectively, and C_{N} is the node capacitance.

Considering this time, the full SRAM cell composed of the two cross-coupled inverters, we obtain the following system of two coupled differential equations:

(9) |

Eq. (8) for the sole and isolated inverter #1 or Eqs. (8) and (9) for the full SRAM cell can be easily solved in the time domain, using a fourth-order Runge-Kutta method [19] with a time step Δt identical to the one used in the RWDD algorithm for charge transport and collection. Using notations introduced in Eq. (8), the incremental results for t_{n+1} = t_{n} + Δt are:

with:

Eqs. (10)–(19) constitute the core model implemented in the code RWDDCPP for CMOS inverter and SRAM cell solving. In Eq. (9), NMOS and PMOS source-to-drain currents are analytically modeled using the EPFL-EKV model 2.6 [20, 21], here implemented in C++ as numerical functions (I_{DN} and I_{DP}) having two arguments (see **Figure 5**). The EKV 2.6 MOSFET model is a predictive (scalable) compact model for the simulation of submicron CMOS technologies. It was built taking into account fundamental physical characteristics of the MOS structure. The model offers a continuous modeling of the different regimes of the transistor operation (weak, moderate and strong inversion), which is mandatory for circuit modeling. The version 2.6 takes into account numerous essential issues for the transistor modeling such as the effects of the doping profile, substrate effects, process-related aspects (oxide thickness, effective channel length and width and junction depth), mobility effects due to vertical and lateral electric fields, the velocity saturation, short-channel effects, etc. Section 5 details several examples of simulations using this code RWDDCPP.

### 4.2. SPICE simulator

For more complex circuit architectures than a single inverter or a SRAM cell, an external SPICE circuit simulator can be otherwise used in the place of the internal subroutine for circuit solving. In this case, the circuit simulator is instantiated in interactive mode by the executable code with a circuit netlist corresponding to the simulated structure. Current source(s) emulating the transient collected current(s) at the different circuit nodes is (are) automatically added to the netlist by the RWDDCPP program. At each time step of the simulation, the magnitude of each current source is updated as a function of the number of charge packets collected by the corresponding biased junction, following Eq. (5). Then, the SPICE simulator computes the new voltages on the circuit nodes and finally, from Eqs. (6) and (7), the width of the space charge region(s) and the distribution of the electrical field can be updated, time step by time step. The whole circuit counteraction to the radiation transient can then be taken into account. This approach has been successfully implemented in this work using both NGSPICE [22] and ELDO [23] simulators.

## 5. Simulation results on various circuits

Four typical simulation cases are presented in this section to validate the proposed model and its computational implementation. All examples have been derived from the same generic 65 nm CMOS technology. In the first example, the RWDD transient response of an isolated transistor subjected to an alpha particle has been compared to TCAD simulation. In the second example, a CMOS inverter has been considered to the effect of the feedback of the node voltage on the modulation of the space charge region related to the impacted junction. The two last examples concern a SRAM cell and a master-slave D flip-flop, with an isolated output, illustrating the internal/external solutions for circuit solving and charge-sharing effects at the level of circuit sensitive nodes.

### 5.1. Model validation—isolated NMOS

This simplest case corresponds to an isolated transistor subjected to single event irradiation. The radiation-induced collected current and charge have been separately computed using the RWDDCPP code and the commercial Synopsys TCAD simulation platform [24]. **Figure 6** (top) shows the 3D structure of a single NMOS transistor implanted in a Pwell region (delimited in depth by a deep N-well) built in Synopsys Sentaurus using geometrical and technological parameters corresponding to the considered 65 nm CMOS technology.

**Figure 6** also shows the electron density 2 ps after the passage of an alpha particle (incident energy of 5 MeV) at the level of the drain of the transistor maintained in the OFF-state. A similar structure has been built and simulated with the RWDDCPP program; the main simulation parameters and options in both approaches have been forced to the same values: drift-diffusion model, constant carrier mobility for electrons (400 cm^{2}/V/s) and holes (200 cm^{2}/V/s), identical minority carrier lifetime (10^{−8} s), etc. **Figure 6** (bottom) shows the transient current and the corresponding charge collected by the drain of the NMOS transistor subjected to a 5 MeV alpha particle striking the center of the drain perpendicularly to the device. A good agreement between results obtained with RWDD and TCAD are shown in this graph. One can observe a few differences between the two curves, especially in the early stages of the transient, mainly due to: (i) the arrival of discrete charge packets at the level of the collecting drain contact in the RWDD approach which induces inherent granularities; and (ii) the characteristic time (2 ps) of the Gaussian time distribution used in the TCAD simulation. Values of the collected charges obtained from TCAD and RWDD are very close at the end of the transient event, with a difference limited to a maximum of 15%, indicating that both modeling approaches have very similar charge collection efficiencies. We confirmed this point by comparing the results given by RWDD and TCAD concerning the space charge region width and the electric field: a very good agreement between the two series of results was obtained and then reiterated for different bias conditions. This shows that both transport and charge collection processes of the RWDD model provide results very similar to those of TCAD simulations without employing any fitting or unphysical parameter. This is a significant advantage of the RWDD approach over other previously developed methods (see [25] for example).

### 5.2. Transient simulation of a CMOS inverter

We explored next the case of the same NMOS transistor but embedded in a CMOS inverter, as shown in **Figure 7**. In the present case, the voltage node V_{2} is not fixed, as was the case for the study of a NMOS transistor standalone; the operation of the second transistor of the CMOS inverter (see **Figure 7**) leads to the variation of the voltage V_{2} during the transient process. These node bias changes may modify the charge collection efficiency at the drain of the NMOS transistor, through the variation of both E and W_{SCR} by a feedback process. In the following, this case is called “with SCR feedback.” By opposition, we name “without SCR feedback” the case when the same (fixed) SCR parameters are maintained during the entire transient simulation. With respect to the charge collection process, this last case is equivalent to that of considering a standalone NMOS transistor. Simulated transients of the CMOS inverter (under the initial conditions V_{1} = 0 and V_{2} = V_{DD}) obtained with the RWDD model when a 5 MeV alpha particle passes across the NMOS drain are shown in **Figure 8**. These results reveal the influence of the circuit feedback on the SCR characteristics, through the time changes of V_{2}. Without the circuit feedback on the SCR, the radiation effect is to turn V_{2} to negative values for a period of time equal to about 30 ps. This obviously evidences a hazardous condition for a cross-coupled inverter in the case of an SRAM cell (somewhat reduced by the gate capacitance in the case where the second inverter is coupled). **Figure 8** illustrates an important issue concerning the inverter operation: the hazardous condition described above persists 2× as long when SCR feedback is activated.

**Figure 8** also shows that the shape of the transient current is essentially modified when the “SCR feedback” is taken into account. The results show a plateau on the collected current curve, in good accordance with TCAD simulations published in [2]. In the case “without feedback,” the transient current shows a classical shape (being composed from a fast drift component and a slower decay) because the SCR of the struck NMOS is maintained constant without feedback. When the “SCR feedback” is considered, the SCR of the struck NMOS can vary when V_{2} changes. As explained by Ferlet-Cavrois et al. [2], when the NMOS is disturbed by a single event, it momentarily biases the on-state load PMOS transistor in a condition for drain current to flow. After a short-duration current peak corresponding to the output capacitance support of the drain voltage at its pre-strike value, the drain voltage collapses, which causes a reduction of both the electric field and the extension of the space charge region; as a consequence, the NMOS single event current becomes governed by the depressed drain voltage (that impacts the collection efficiency of the junction) and the compensating PMOS transistor drive current. The result of this dynamic interaction of the node voltage and the PMOS transient current is a characteristic current equilibrium or “plateau,” as shown in **Figure 8**. The amplitude of this current plateau is linked to the PMOS drive, and the duration of the plateau is synchronous with the depressed drain voltage. As soon as most of the deposited charge flows out of the struck transistor, the current flow cannot be maintained, the equilibrium conditions relax, the drain voltage recovers and the current pulse again decreases toward zero [2, 26]. This second example illustrates the importance of device coupling effects at circuit-level and its consequences on the charge collection dynamics in the impacted device. The soft-error rate consequence of such effects has to be investigated in future works, notably in SRAM cells composed of two cross-coupled inverters for which a very similar behavior as highlighted in **Figure 8** is expected at the level of the impacted OFF-state transistor.

### 5.3. Alpha-particle SEU in SRAM

In this third example, the single-event upset alpha particle cross-section of a 65 nm SRAM cell has been simulated using the RWDD approach and compared with experiments. The alpha cross-section experimental measurement was performed using an Americium 241 source (3.7 MBq), at room temperature and under nominal voltage 1.2 V), following all the recommendations of the JESD98A test standard [27]. Since the alpha source has a high activity, several thousand errors are measured and cross-section uncertainties are very low.

For the RWDD simulation, several additional assumptions have been introduced to estimate the electrical response of the 65 nm bit cell: (i) all of the electron-hole pairs generated inside or below the deep N-well (see the structure of the transistor in **Figure 6** top) are not taken into account; (ii) minorities charge carriers that reach well limits are considered as recombined and are eliminated from the simulation. **Figure 9** (top) shows the 65 nm bit cell composed of a centered Nwell and two Pwells. The four drains of the transistors constituting the two inverters have been placed in structure and connected to the corresponding nodes in the circuit netlist. The total simulation area is larger than the bit cell area to take into account an alpha particle with a high tilt angle, which strikes around the bit cell. The maximum tilt angle has been calculated using the range of the alpha particle and the thickness of the back-end layer. We computed a simulated cross-section by performing simulations for several thousand impacts of alpha particles, considering random impact locations and particles with random directions; the number of upsets was divided by the simulated fluency.

The measured alpha particle cross-section of the 65 nm bit cell is compared with their simulated counterparts in **Figure 9** (bottom). Three different voltages from 1.08 to 1.32 V have been considered in simulation. This figure illustrates the nice agreement between the results obtained with the RWDD model and experimental measurements. **Figure 9** also demonstrates that this model is able to take into account the supply voltage dependence for the evaluation of the bit cell cross-section. Unfortunately, the comparison between experiment and simulation was not possible for this bit cell at both 1.1 and 1.3 V because experimental data were not available.

### 5.4. Flip-flop circuit

The final step was to confirm the ability of the RWDD model to simulate more complex circuit architectures and considering a large LET range. For this purpose, we simulated the heavy ion cross-section of a 65 nm flip-flop and we compared the simulation results with experimental measurements. To model the geometrical structure of the flip-flop, we took into account simulation hypothesis identical to those considered for the SRAM. Heavy ion measurements have been performed on a 65 nm test-chip, at nominal supply voltage and room temperature. The measurements have been performed at the RADEF facility, Jyväskylä University (Finland), with a fluency of 10^{6} heavy ions per square centimeter for five different ion LETs.

The functional schematic of the considered flip-flop (classical, unhardened master slave D-flip-flops) is shown in **Figure 10** (top). Data are transmitted between the master and slave latch using a pass gate. In this study, we include in the modeled structure all drains of all the transistors present in the flip-flop, unlike to our previous simulation works that only considered the bistable transistors [6, 28]. In this way, all of the possible upset mechanisms can be taken into account. This is specifically the case when an ionizing particle strikes a clock network transistor that can involve the latching of new data, these mechanisms being particularly important for high LET ions.

**Figure 10** (bottom) shows the measured and simulated heavy ion cross-section of the flip-flop. This graph shows the excellent correlation of the simulated cross-section with the measurements for low LET ions. For higher LET impacts, the non-modeled multi-cell upsets induces a difference between measured and simulated data. Several reasons can explain these differences:

- The main reason for the lack of accuracy is the contribution of non-modeled multi-cell upset (MCU) in the measured data: indeed, flip-flop cross-sections have been measured in a shift register structure where flip-flops are next to each other. Thus, a single ionizing particle can upset several flip-flops. Experimentally, this phenomenon has been verified (some MCUs have been recorded) but it was not possible to precisely quantify it since the number of stages in the shifter structure is not constant between the physical neighbor’s cells, and the test has been performed dynamically.

- A second reason that can explain the differences in the cross-section at high LET energies concerns modeling of the charge carrier reaching the good limits and more specifically the approximation used to calculate this charge. Charge packets that reach the good interface are considered in the charge transport model as recombined and are then removed from the simulation. Another limitation concerns the impact location of the ionizing particle: the current implementation of our model is not totally adequate for an impact close to the good interface. Owing to the punctual character of the track profile, all charge packets are generated on the impacted well since, if taking into account the particle track radius, a non-negligible part of these packets should be considered in the adjacent well(s).

The SEU sensitive areas for each simulated LET can be identified on the flip-flop layout using the RWDD model, as in the case for other Monte Carlo simulation methods [29, 30]. The map of the simulated impacts inducing an upset for each ion LET is illustrated in **Figure 11**. This map has been generated on the slave latch of the flip-flop, with a high logic state stored in the latch. This figure indicates that, as the LET of the ionizing particles increases, the sensitive transistors are more abundant. The explanation is that upsets are uniquely produced on bistable transistors for ions with low LET, while for higher ion LET, other transistors (clock and pass gate transistors) become sensitive to SEU, in addition to those impacted at low LET.

Finally, in order to evaluate the capability of the RWDD model to simulate charge-sharing mechanisms, a dedicated study has been performed. The SEU impact map of **Figure 11** shows that the top part of the bottom right drain of the flip-flop (labeled node “SL”) is not sensitive to SEU. A detailed study of the collected currents in this area has then been conducted. For this purpose, a specific ion impact location (shown in **Figure 11**, labeled “detailed impact” in the top left corner of node SL) has been fixed and the collected current waveform following the ion impact at this location has been systematically analyzed, as well as the voltage waveform of the two bistable nodes. As expected, this node SL collects most of the charge as it corresponds to the stricken drain. However, the fine analysis (not shown) of the different transient currents indicates that the nodes SLN, N1 and N2 collect later, by a diffusion mechanism, a small quantity of charge, which is high enough to keep the voltage of the SLN node at the low logic state and to prevent a latch upset. In order to validate this affirmation, we removed from the circuit netlist in a separate simulation batch the drains corresponding to nodes SLN, N1 and N2. In this case, since no charge is collected on the NMOS of the SLN node to counter-balance the feedback loop of the latch, the electrical potential of the SLN node increases and an upset occurs.

To conclude, this particular charge-sharing mechanism in the area of the top left corner of collecting node SL has been observed whatever the value of the LET in the range 1.8–60 MeV cm^{2}/mg, demonstrating the capability of the RWDD approach to treat such complex circuit response to single events.

## 6. Conclusion

A new computational method for the simulation of single event effects in integrated circuits has been presented. This approach is a Monte Carlo method based on a random-walk drift-diffusion algorithm that transports the radiation-induced charge, segmented into discrete charge packets, in the semiconductor regions of a given circuit architecture. Carrier diffusion is very well reproduced with the random-walk algorithm while the carrier “drift” component of the model perfectly captures the effects of the electric field developed in the space-charge region of the reversely biased collecting contact(s). The model has been fully derived in C++ using advanced structures to model device/circuit geometry, particle track and charge transport, collection, recombination and extraction. This approach has been dynamically coupled with an internal subroutine or an external circuit simulator to take into account spatial and temporal variations of the electric field in the vicinity of the collecting structure(s). Thus, complex architectures, such as flip-flops, can be easily modeled and charge-sharing mechanisms are accurately simulated.

This chapter mainly focused on the model implementation and the way to solve the circuit response in the time domain, taking into account the circuit feedback on the charge collection process. Four simulation test cases have been explored and compared to radiation experiments or TCAD simulations in order to validate the proposed model. These test cases show good quantitative agreement between measurements and simulated data over a large range of LETs up to 60 MeV.cm^{2}/mg and structure complexity. This first implementation remains therefore not self-consistent with the electrostatic potential (in other words, Poisson’s equation is not solved during the transient computation) and does not take into account the possible interactions between charge packets. The model ability to accurately reproduce a regime of high carrier injection is therefore uncertain. The limitations of the method will be more quantitatively explored in a forthcoming dedicated study, as well as possible improvements in terms of self-consistency with the electrostatic problem.