This example describes the simulation of a silicon diode. We will look at the carrier profiles and lifetimes as a function of bias/doping as well as I-V characteristics in both reverse and forward bias under several bulk recombination processes. In particular, we will inspect the effect of band-to-band tunneling and impact ionization on the current.
Simulation setup
Necessary project files and a script file are provided to assist with the application example. The project files contain the material properties, geometry, and simulation region required to run the example.
Open the silicon_diode.ldev file in CHARGE. Notice that the diode is specified as follows:
- The thickness of the active silicon layer is 4um.
- Aluminum is used for both the cathode and anode contacts.
- By default, trap assisted (Shockley-Read-Hall) recombination, Auger recombination and radiative recombination are turned on.
- To define the space-charge regions in the silicon, first the background doping concentration of the silicon is set to represent a p-type epitaxial layer with a concentration of 1x1018 cm-3. This is accomplished by defining a region of constant doping that encompasses the entire geometry. Next, a diffusion doping region is used to specify the n-dopant concentration under the cathode with a peak doping concentration of 1x1020 cm-3 at the surface.
Doping
The doping concentration can be visualized as shown below. Click on the "mesh" button to mesh the simulation region. Once the meshing is complete, right click on CHARGE object and visualize "grid" where you can get an image plot of N (the net dopant concentration). Alternatively, you can use the script below to interpolate the finite element mesh results onto a linear mesh for a 1D plot:
v_index=1; x = -0.025e-6; z = linspace(getnamed('silicon','z min'),getnamed('cathode','z min'),401); doping = getresult('CHARGE','grid'); vtx = [doping.x, doping.z]; N = abs(doping.N); N_src = interptri(doping.elements,vtx,N,x,z,0); plot(z*1e6,N_src,'depth (um)','doping (1/cm3)','', 'logplot');
Material Properties
In this example two materials, silicon and aluminum are used. Expand the Material Group in the Objects Tree. Here, you can view the materials used in the simulation. Two types of materials are used in this example:
- Conductor (Aluminum for the anode and cathode contacts)
- Semiconductor (Silicon for the bulk silicon)
The conductor is defined by its work function. The semiconductor is defined using the more complex semiconductor models. In particular the "Recombination" tab of this material is of interest in this example.
For all parts of this example, the bulk recombination models, trap-assisted (SRH), radiative, and Auger recombination are enabled.
Trap assisted tunneling can be enabled using the field drop down menu under "Carrier Lifetime Correction Models" section (see image to the right). When enabled, two models, Hurkx and Schenk are available. For the Hurkx model used here the tunneling mass coefficient mt is 0.25.
Band to band tunneling can be enabled using the corresponding tick-box in the "model" list. When enabled, two models, Hurkx and Schenk are available. For the Hurkx model used here, the coefficients are as follows:
- A = 4e14 cm-3s-1
- B = 1.9e7 V/cm
- Eta =2.5
Similarly, for "Impact Ionization", the user can specify a model for impact ionization of electrons and holes using the Selberherr model. For this model the default coefficients for silicon are as follows:
Coefficient |
Electrons |
Holes |
---|---|---|
Ecrit (V/cm) |
1.23x106 |
2.03x106 |
Alpha (1/cm) |
7.03x105 |
1.58x106 |
Beta |
1 |
1 |
For descriptions of the recombination processes as well as the tunneling models, please see Lumerical's knowledgebase entry on Semiconductors.
Different recombination models in the semiconductor material properties
In the silicon_diode.ldev project, check the Si material properties to make sure that all three recombination processes (radiative, Auger, and trap assisted) are enabled. Also verify that trap assisted tunneling, band to band tunneling, and impact ionization are disabled. Set the contact biases in the Boundary Conditions group according to the table below. You can edit the properties of the anode and cathode boundary conditions (electrical contacts) by selecting them and clicking on the "Edit" button. By sweeping the cathode bias over a range of values, we are reverse biasing the diode.
contact |
Voltage (volts) |
---|---|
cathode |
0-2 in 6 steps |
anode |
0 |
Run and results
Run the simulation. Once the simulation is run, you can right click on the CHARGE object, and visualize the cathode dataset. In order to get only the relevant I plot, you can remove all the attributes except for I. You should get the plot below:
Reverese bias current of the silicon diode
You can also visualize other quantities such as variation of trap assisted recombination with voltage or the doping dependance of carrier lifetimes. In the material properties, for silicon, note that the carrier lifetimes are set to the same value for both electrons and holes and in the correction models become a function of doping concentrations. A script similar to the one above can be used to interpolate the carrier lifetime and recombination rate onto a linear plot. Run the script file silicon_diode_linear_plots.lsf to plot the electron lifetime (left, below). Note that the electron (and hole) lifetime follows the same pattern as the doping profile. The script also plots the trap assisted recombination rate for 4 different bias voltages. At equilibrium, recombination is negligible but becomes larger as the bias is increased. The curve has been zoomed on to the space-charge region. The width of the peak increases with bias which also indicates the increased width of the space charge region.
Doping dependence of carrier lifetime
Trap assisted recombination rate at different reverse bias voltages
Trap assisted recombination
Next, we want to model the device in forward bias. Set the contact biases according to the table below.
contact |
Voltage (volts) |
---|---|
cathode |
0 |
anode |
0-0.4 in 5 steps |
Under the sweeps and optimizations tab, a sweep task has been set up to sweep over the temperature of the simulation from 300 to 360 K in three points. The anode current is extracted. Right click on the sweep task and click run to run the sweep. This will run three files with 5 bias points in each. Once the sweep is done, run the following lines of script to extract the current from the sweep and plot the results as shown left below.
T = getsweepresult('sweep','T'); Ia = getsweepresult('sweep','Ia'); I = pinch(Ia.I); V = pinch(Ia.V_anode); plot(V(2:5),1e6*I(2:5,1:3),'voltage (V)', 'current (uA)','','logplot'); legend('T=300K','T=330K','T=360K');
I-V characteristics at different temperature with no trap assisted tunneling
I-V characteristics at different temperature including trap assisted tunneling
Now, update the silicon material model to account for trap assisted tunneling and repeat the sweep again. To do this, from the field drop down menu under trap assisted recombination, choose the Hurkx model. Verify that the effective tunneling mass coefficient mt is 0.25. This creates a temperature ( and doping) dependent carrier lifetime correction model. Run the sweep, and using the script above, the plot on the right (above) will be generated.
Band to band tunneling and impact ionization
Next, we will put the diode back in reverse bias mode. This time, we want to observe the effect of band to band tunneling and impact ionization. Use the silicon_diode_bbt.ldev file for this part of the example.
NOTE: field-dependent breakdown models such as band to band tunneling and impact ionization are exponentially sensitive to variations in the electrostatic field and potential, which can lead to convergence failure. To improve the convergence of the simulation, we recommend increasing the iteration limit and enabling gradient mixing (both under the advanced settings tab of the CHARGE solver). Additionally, it is important to accurately resolve the electric fields in the critical regions, and mesh overrides may be necessary to ensure that the mesh is sufficiently refined; edge lengths of 4 nm or less are often required. |
To accurately resolve the electric fields and ensure the convergence of the simulation, we have added mesh override regions at the junction and set the iteration limit under the advanced tab of the CHARGE solver region to 240. To improve convergence, we also set the gradient mixing setting to "fast". Because we are increasing the mesh density, we will reduce the x span to improve simulation performance (make it faster).
Set the contact biases (boundary conditions) according to the table below. The number of voltage steps to sweep from 0 to 3 volts is also increased to improve the convergence of the simulation.
contact |
Voltage (volts) |
---|---|
cathode |
0-3 v in 61 steps |
anode |
0 |
In this example, we will run three different simulations with the same bias settings: default material settings with no band to band tunneling or impact ionization, silicon with band to band tunneling only, and silicon with both band to band tunneling and impact ionization.
First, run the simulation with no band to band tunneling or impact ionization and for band to band tunneling only. Once the simulations are done, use the following script to plot the current in log scale. Note that we have set the plot range in y axis from 1e-12 to 1e-8 to more easily view the results. The two plots below are for no tunneling (left) and band to band tunneling (right) only.
anode = getresult('CHARGE','anode'); I = abs(pinch(anode.I)); V = anode.V_cathode; plot(V,I, 'Vanode (v)', 'Ianode (A)', '','logplot'); setplot('y min',1e-12); setplot('y max',1e-8);
Reverse bias current with no tunneling
Reverse bias current with band to band tunneling
As you run the impact ionization case, you will notice that the solver fails to converge at a bias of 2 volts, which indicates the impact ionization break down bias for the device. At this point the current value increases up by several orders of magnitude for small changes in the voltage (~1mV). Choose "quit and save" to save the results and the script below to plot the current for the range 0 to 2 volts. Alternatively, re-run the simulation with the cathode voltage going from 0 to 2 volt in 41 steps and use the script above to plot the current.
anode = getresult('CHARGE','anode'); bbtii = abs(pinch(anode.I)); V1 = linspace(0,2,41); lastgood = 41; plot(V1,bbtii(1:lastgood),'Vanode (v)', 'Ianode (A)', '','logplot'); setplot('y min',1e-12); setplot('y max',1e-8);
Note: The critical field parameter for the impact ionization model was intentionally reduced to artificially lower the breakdown voltage in the silicon_diode_bbt.ldev file so that we can observe the phenomenon faster. The default parameters are listed in the table at the beginning of this example. |
Here we have plotted all three I-V curves on the same plot for comparison purposes. You can plot the results from the three simulations together using the holdon, holdoff commands.
Effect of different recombination mechanism on reverse bias current