ViennaSHE 1.3.0
Free open-source semiconductor device simulator using spherical harmonics expansions techniques.
nin-diode-1d.cpp
Go to the documentation of this file.
1/* ============================================================================
2 Copyright (c) 2011-2022, Institute for Microelectronics,
3 Institute for Analysis and Scientific Computing,
4 TU Wien.
5
6 -----------------
7 ViennaSHE - The Vienna Spherical Harmonics Expansion Boltzmann Solver
8 -----------------
9
10 http://viennashe.sourceforge.net/
11
12 License: MIT (X11), see file LICENSE in the base directory
13=============================================================================== */
14
15#if defined(_MSC_VER)
16 // Disable name truncation warning obtained in Visual Studio
17 #pragma warning(disable:4503)
18#endif
19
20#include <iostream>
21#include <cstdlib>
22#include <vector>
23
24// ViennaSHE includes:
25#include "viennashe/core.hpp"
26
27// ViennaGrid mesh configurations:
28#include "viennagrid/config/default_configs.hpp"
29
30
58template <typename DeviceType>
59void init_device(DeviceType & device)
60{
61 typedef typename DeviceType::segment_type SegmentType;
62
64 SegmentType const & contact_left = device.segment(0);
65 SegmentType const & n_left = device.segment(1);
66 SegmentType const & i_center = device.segment(2);
67 SegmentType const & n_right = device.segment(3);
68 SegmentType const & contact_right = device.segment(4);
69
71 device.set_material(viennashe::materials::si(), n_left);
72 device.set_material(viennashe::materials::si(), i_center);
73 device.set_material(viennashe::materials::si(), n_right);
74
75 // Contacts are always conductors/metals (for now)
76 device.set_material(viennashe::materials::metal(), contact_left);
77 device.set_material(viennashe::materials::metal(), contact_right);
78
84 std::cout << "* init_device(): Setting doping..." << std::endl;
85 device.set_doping_n(1e25, n_left);
86 device.set_doping_p(1e6, n_left);
87
88 device.set_doping_n(1e20, i_center);
89 device.set_doping_p(1e11, i_center);
90
91 device.set_doping_n(1e25, n_right);
92 device.set_doping_p(1e6, n_right);
93
95 device.set_contact_potential(0.0, contact_left);
96 device.set_contact_potential(0.5, contact_right);
97}
98
107int main()
108{
114
115 std::cout << viennashe::preamble() << std::endl;
116
118 std::cout << "* main(): Creating device..." << std::endl;
119 DeviceType device;
120
128
137 // Left contact
138 generator_params.add_segment(0.0, 1e-6, 10); //start at x=0, length 1e-6, 10 points
139 // n_left
140 generator_params.add_segment(1e-6, 1e-6, 10);
141 // i_center
142 generator_params.add_segment(2e-6, 1e-6, 10);
143 // n_right
144 generator_params.add_segment(3e-6, 1e-6, 10);
145 // Right contact
146 generator_params.add_segment(4e-6, 1e-6, 10);
147
149 device.generate_mesh(generator_params);
150
154 std::cout << "* main(): Initializing device..." << std::endl;
155 init_device(device);
156
157
164 std::cout << "* main(): Creating DD simulator..." << std::endl;
165
166
174 viennashe::config dd_cfg;
175
176 // enable electrons and holes, use the continuity equation for each of them:
177 dd_cfg.with_holes(true);
178 dd_cfg.with_electrons(true);
179
182
183 // The default nonlinear solver is a Gummel iteration. Here we use Newton's method with at most 50 iterations:
185 dd_cfg.nonlinear_solver().max_iters(50);
186
193 viennashe::simulator<DeviceType> dd_simulator(device, dd_cfg);
194
195 std::cout << "* main(): Launching simulator..." << std::endl;
196 dd_simulator.run();
197
203 viennashe::io::write_quantities_to_VTK_file(dd_simulator, "nin1d_dd_quan");
204
205
211 std::cout << "* main(): Setting up SHE..." << std::endl;
212
219 viennashe::config config;
220
221 // Config a first-order SHE for electrons and holes
222 config.max_expansion_order(1);
223 config.with_electrons(true);
224 config.with_holes(true);
227
228 // Scattering mechanisms for SHE: We only consider phonon scattering:
229 config.scattering().acoustic_phonon().enabled(true);
230 config.scattering().optical_phonon().enabled(true);
231 config.scattering().ionized_impurity().enabled(false);
232 config.scattering().impact_ionization().enabled(false);
233 config.scattering().electron_electron(false);
234
235 // nonlinear solver: At most 30 Gummel iterations with a damping of 0.3:
236 config.nonlinear_solver().max_iters(30);
237 config.nonlinear_solver().damping(0.4);
238
239 // Configure SHE: Minimum energy range is 0.5 eV
241
242 // energy spacing of 31meV
244
245
253 std::cout << "* main(): Computing SHE..." << std::endl;
254 viennashe::simulator<DeviceType> she_simulator(device, config);
255
256 // Use DD solution as an initial guess
257 she_simulator.set_initial_guess(viennashe::quantity::potential(), dd_simulator.potential());
258 she_simulator.set_initial_guess(viennashe::quantity::electron_density(), dd_simulator.electron_density());
259 she_simulator.set_initial_guess(viennashe::quantity::hole_density(), dd_simulator.hole_density());
260
261 she_simulator.run();
262
268 std::cout << "* main(): Writing SHE result..." << std::endl;
269
271 she_simulator.config(),
272 she_simulator.quantities().electron_distribution_function(),
273 "nin1d_edf");
274
277 viennashe::io::write_quantities_to_VTK_file(she_simulator, "nin1d_she_quan");
278
279
283
284 // Helper type returning the generalized energy distribution function:
286
287 // Write the generalized EDF to a plain text file
288 gedfwriter(device,
289 viennashe::util::any_filter(), // write the result for all cells
290 EDFWrapperType(she_simulator.config(), she_simulator.quantities().electron_distribution_function()),
291 "nin1d_edf.dat"); // filename
292
293 // Write the electron density to a plain text file
295 gpwriter(device,
296 viennashe::util::any_filter(), // write the result for all cells
297 she_simulator.electron_density(),
298 "nin1d_n.dat"); // filename
299
301 std::cout << "* main(): Results can now be viewed with your favorite VTK viewer (e.g. ParaView)." << std::endl;
302 std::cout << "* main(): Don't forget to scale the z-axis by about a factor of 1e12 when examining the distribution function." << std::endl;
303 std::cout << std::endl;
304 std::cout << "*********************************************************" << std::endl;
305 std::cout << "* ViennaSHE finished successfully *" << std::endl;
306 std::cout << "*********************************************************" << std::endl;
307
308 return EXIT_SUCCESS;
309}
310
The main SHE configuration class. To be adjusted by the user for his/her needs.
Definition: config.hpp:124
double min_kinetic_energy_range(viennashe::carrier_type_id ctype) const
Returns the minimum kinetic energy range for the selected carrier.
Definition: config.hpp:388
long max_expansion_order() const
Returns the current maximum expansion order.
Definition: config.hpp:369
bool with_holes() const
Returns true if holes are considered in the simulation.
Definition: config.hpp:234
nonlinear_solver_config_type & nonlinear_solver()
Returns the configuration object for the nonlinear solver.
Definition: config.hpp:495
void set_electron_equation(equation_id equ_id)
Definition: config.hpp:231
double energy_spacing() const
Returns the uniform energy spacing of discrete energies.
Definition: config.hpp:460
void set_hole_equation(equation_id equ_id)
Definition: config.hpp:239
bool with_electrons() const
Returns true if electrons are considered in the simulation.
Definition: config.hpp:226
viennashe::she::scatter_config & scattering()
Returns the configuration object for scattering.
Definition: config.hpp:475
Defines the physical properties of a device, e.g. doping. This is the implementation for 2d and highe...
Definition: device.hpp:818
A convenience wrapper for accessing the generalized energy distribution function (f * Z,...
ionized_impurity_scattering_parameters const & ionized_impurity() const
Returns the parameters for ionized impurity scattering. Const-version.
Definition: config.hpp:377
impact_ionization_scattering_parameters const & impact_ionization() const
Returns the parameters for impact ionization scattering. Const-version.
Definition: config.hpp:382
acoustic_phonon_scattering_parameters const & acoustic_phonon() const
Returns the parameters for acoustic phonon scattering. Const-version.
Definition: config.hpp:367
optical_phonon_scattering_parameters const & optical_phonon() const
Returns the parameters for optical phonon scattering. Const-version.
Definition: config.hpp:372
bool electron_electron() const
Returns true if electron-electron scattering is activated.
Definition: config.hpp:387
Class for self-consistent SHE simulations.
Definition: simulator.hpp:675
std::size_t max_iters() const
Returns the maximum number of nonlinear solver iterations.
Definition: config.hpp:317
double damping() const
Returns the damping for the nonlinear solver.
Definition: config.hpp:328
void set(long solver_id)
Sets a new linear solver. Provide IDs defined in.
Definition: config.hpp:270
Configuration class for the simple mesh generator.
void add_segment(double start_x, double start_y, double len_x, double len_y, unsigned long points_x, unsigned long points_y)
Convenience header, which includes all core functionality available in ViennaSHE.
void init_device(DeviceType &device, double Vg_init, double Nd, double Na)
Initalizes the MOS 1D device for testing.
Definition: mos1d_dg.hpp:94
void write_quantities_to_VTK_file(viennashe::simulator< DeviceType > const &simulator_obj, std::string filename, bool include_debug_information=false)
Generic interface function for writing simulated quantities to a VTK file.
std::string electron_density()
@ EQUATION_SHE
Definition: forwards.h:118
@ EQUATION_CONTINUITY
Definition: forwards.h:117
std::string preamble()
Prints the ViennaSHE preamble (header). Used in all the examples as well as the standalone-applicatio...
Definition: version.hpp:63
int main()
Definition: resistor1d-c.c:108
Writes the energy distribution function to a file which can be processed by Gnuplot....
Writes quantities to a file which can be processed by Gnuplot. Works for 1d, 2d and 3d only.
A class referring to any metal contact.
Definition: all.hpp:44
A class referring to silicon and providing certain material parameters Note that this is the default ...
Definition: all.hpp:50
static const double q
Elementary charge.
Definition: constants.hpp:44
@ newton_nonlinear_solver
Gummel iteration.
Definition: config.hpp:244
A trivial filter, all objects are accepted.
Definition: filter.hpp:41