ViennaSHE 1.3.0
Free open-source semiconductor device simulator using spherical harmonics expansions techniques.
mosfet.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
67template <typename DeviceType>
68void init_device(DeviceType & device)
69{
70 typedef typename DeviceType::segment_type SegmentType;
71
73 SegmentType const & gate_contact = device.segment(1);
74 SegmentType const & source_contact = device.segment(2);
75 SegmentType const & gate_oxide = device.segment(3);
76 SegmentType const & drain_contact = device.segment(4);
77 SegmentType const & source = device.segment(5);
78 SegmentType const & drain = device.segment(6);
79 SegmentType const & body = device.segment(7);
80 SegmentType const & body_contact = device.segment(8);
81
83 std::cout << "* init_device(): Setting materials..." << std::endl;
84 device.set_material(viennashe::materials::metal(), gate_contact);
85 device.set_material(viennashe::materials::metal(), source_contact);
86 device.set_material(viennashe::materials::metal(), drain_contact);
87 device.set_material(viennashe::materials::metal(), body_contact);
88
89 device.set_material(viennashe::materials::hfo2(), gate_oxide);
90
91 device.set_material(viennashe::materials::si(), source);
92 device.set_material(viennashe::materials::si(), drain);
93 device.set_material(viennashe::materials::si(), body);
94
101 std::cout << "* init_device(): Setting doping..." << std::endl;
102 device.set_doping_n(1e24, source);
103 device.set_doping_p(1e8, source);
104
105 device.set_doping_n(1e24, drain);
106 device.set_doping_p(1e8, drain);
107
108 device.set_doping_n(1e17, body);
109 device.set_doping_p(1e15, body);
110
115 device.set_contact_potential(0.8, gate_contact);
116 device.set_contact_potential(0.0, source_contact);
117 device.set_contact_potential(1.0, drain_contact);
118 device.set_contact_potential(0.0, body_contact);
119
120}
121
122
131int main()
132{
138 typedef DeviceType::segment_type SegmentType;
139
140 std::cout << viennashe::preamble() << std::endl;
141
148 std::cout << "* main(): Creating and scaling device..." << std::endl;
149 DeviceType device;
150 device.load_mesh("../examples/data/mosfet840.mesh");
151 device.scale(1e-9);
152
153
157 std::cout << "* main(): Initializing device..." << std::endl;
158 init_device(device);
159
160
169 std::cout << "* main(): Creating DD simulator..." << std::endl;
170
178 viennashe::config dd_cfg;
179
180 // enable electrons and holes and specify that for each of them a continuity equation should be solved:
181 dd_cfg.with_electrons(true);
182 dd_cfg.with_holes(true);
185
186 // Nonlinear solver parameters: 200 Gummel iterations with rather high damping
187 dd_cfg.nonlinear_solver().max_iters(200);
188 dd_cfg.nonlinear_solver().damping(0.125);
189
190
197 viennashe::simulator<DeviceType> dd_simulator(device, dd_cfg);
198 std::cout << "* main(): Launching DD simulator..." << std::endl;
199 dd_simulator.run();
200
208 viennashe::io::write_quantities_to_VTK_file(dd_simulator, "mosfet_dd_quan");
209
214 SegmentType const & drain_contact = device.segment(4);
215 SegmentType const & body = device.segment(6);
216
217 std::cout << "* main(): Drain electron current Id_e = " << viennashe::get_terminal_current(device, viennashe::ELECTRON_TYPE_ID,
218 dd_simulator.potential(), dd_simulator.electron_density(),
220 body, drain_contact ) * 1e-6 << std::endl;
221 std::cout << "* main(): Drain hole current Id_h = " << viennashe::get_terminal_current(device, viennashe::HOLE_TYPE_ID,
222 dd_simulator.potential(), dd_simulator.hole_density(),
224 body, drain_contact ) * 1e-6 << std::endl;
225
226
227
241 std::cout << "* main(): Setting up first-order SHE (semi-self-consistent using 40 Gummel iterations)..." << std::endl;
242 viennashe::config config;
243
244 // Set the expansion order to 1 for SHE
245 config.max_expansion_order(1);
246
247 // Use both carrier types
248 config.with_electrons(true); config.with_holes(true);
249
250 // Configure equation
251 config.set_electron_equation(viennashe::EQUATION_SHE); // SHE for electrons
253
254 // Energy spacing of 31 meV
255 config.energy_spacing(31.0 * viennashe::physics::constants::q / 1000.0);
256
257 // Set the scattering mechanisms
258 config.scattering().acoustic_phonon().enabled(true);
259 config.scattering().optical_phonon().enabled(true);
260 config.scattering().ionized_impurity().enabled(false);
261
262 // The linear solver should run for at most 2000 iterations:
263 config.linear_solver().max_iters(2000);
264
265 // Configure the nonlinear solver to 40 Gummel iterations with a damping parameter 0.4
266 config.nonlinear_solver().max_iters(40);
267 config.nonlinear_solver().damping(0.4);
268
276 std::cout << "* main(): Computing first-order SHE..." << std::endl;
277 viennashe::simulator<DeviceType> she_simulator(device, config);
278
279 // Set the previous DD solution as an initial guess
280 she_simulator.set_initial_guess(viennashe::quantity::potential(), dd_simulator.potential());
281 she_simulator.set_initial_guess(viennashe::quantity::electron_density(), dd_simulator.electron_density());
282 she_simulator.set_initial_guess(viennashe::quantity::hole_density(), dd_simulator.hole_density());
283
284 // Trigger the actual simulation:
285 she_simulator.run();
286
292 std::cout << "* main(): Writing energy distribution function from first-order SHE result..." << std::endl;
294 she_simulator.config(),
295 she_simulator.quantities().electron_distribution_function(),
296 "mosfet_she_edf");
297
299 viennashe::io::write_quantity_to_VTK_file(she_simulator.potential(), device, "mosfet_she_potential");
300 viennashe::io::write_quantity_to_VTK_file(she_simulator.electron_density(), device, "mosfet_she_electrons");
301
304 viennashe::io::write_quantities_to_VTK_file(she_simulator, "mosfet_she_quan");
305
310 std::cout << "* main(): Drain hole current Id_e = " << viennashe::get_terminal_current(
311 device, config, she_simulator.quantities().electron_distribution_function(), body, drain_contact ) * 1e-6
312 << std::endl;
313// std::cout << "Id_h = " << viennashe::get_terminal_current(
314// device, config, she_simulator.quantities().hole_distribution_function(), body, drain_contact ) * 1e-6
315// << std::endl;
316
318 std::cout << "* main(): Results can now be viewed with your favorite VTK viewer (e.g. ParaView)." << std::endl;
319 std::cout << "* main(): Don't forget to scale the z-axis by about a factor of 1e12 when examining the distribution function." << std::endl;
320 std::cout << std::endl;
321 std::cout << "*********************************************************" << std::endl;
322 std::cout << "* ViennaSHE finished successfully *" << std::endl;
323 std::cout << "*********************************************************" << std::endl;
324
325 return EXIT_SUCCESS;
326}
327
The main SHE configuration class. To be adjusted by the user for his/her needs.
Definition: config.hpp:124
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
linear_solver_config_type & linear_solver()
Returns the configuration object for the linear solver.
Definition: config.hpp:485
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
ionized_impurity_scattering_parameters const & ionized_impurity() const
Returns the parameters for ionized impurity scattering. Const-version.
Definition: config.hpp:377
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
Class for self-consistent SHE simulations.
Definition: simulator.hpp:675
std::size_t max_iters() const
Returns the maximum number of iterative linear solver iterations.
Definition: config.hpp:155
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
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.
void write_quantity_to_VTK_file(QuantityType const &quantity, DeviceType const &device, std::string filename, std::string name_in_file="viennashe_quantity")
Generic interface function for writing a quantity to a VTK file. Automatically dispatches between ver...
viennashe::models::dd::mobility< DeviceType > create_constant_mobility_model(DeviceType const &device, double mu)
Returns a mobility model, which always yields the same mobility.
Definition: mobility.hpp:61
std::string electron_density()
double get_terminal_current(DeviceT const &device, CurrentAccessorT const &current_accessor, SegmentT const &semiconductor, SegmentT const &terminal)
Returns the terminal current for a number of given vertices. Considers carrier flux and displacement ...
@ HOLE_TYPE_ID
Definition: forwards.h:188
@ ELECTRON_TYPE_ID
Definition: forwards.h:187
@ 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
A class referring to hafnium dioxide.
Definition: all.hpp:119
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