ViennaSHE 1.3.0
Free open-source semiconductor device simulator using spherical harmonics expansions techniques.
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
carrier_density.hpp
Go to the documentation of this file.
1#ifndef VIENNASHE_SHE_CARRIER_DENSITY_HPP
2#define VIENNASHE_SHE_CARRIER_DENSITY_HPP
3
4/* ============================================================================
5 Copyright (c) 2011-2022, Institute for Microelectronics,
6 Institute for Analysis and Scientific Computing,
7 TU Wien.
8
9 -----------------
10 ViennaSHE - The Vienna Spherical Harmonics Expansion Boltzmann Solver
11 -----------------
12
13 http://viennashe.sourceforge.net/
14
15 License: MIT (X11), see file LICENSE in the base directory
16=============================================================================== */
17
18
19#include <math.h>
20#include <fstream>
21#include <iostream>
22
23#include "viennagrid/mesh/mesh.hpp"
24
25#include "viennashe/forwards.h"
26#include "viennashe/config.hpp"
32
37namespace viennashe
38{
39 namespace she
40 {
41
42 namespace detail
43 {
45 template <typename SHEQuantity>
47 {
49
50 public:
51
52 typedef double value_type;
53
55 SHEQuantity const & quan,
56 double energy_start = 0.0,
57 double energy_end = 1.0) : conf_(conf), quan_(quan), energy_start_(energy_start), energy_end_(energy_end) {}
58
60 : quan_(o.quan_), conf_(o.conf_), energy_start_(o.energy_start_), energy_end_(o.energy_end_)
61 { }
62
64 : conf_(o.conf_), quan_(quan), energy_start_(o.energy_start_), energy_end_(o.energy_end_)
65 { }
66
67 template <typename ElementType>
68 value_type operator()(ElementType const & elem) const
69 {
70 value_type density = 0;
71
72 typename viennashe::config::dispersion_relation_type dispersion = conf_.dispersion_relation(quan_.get_carrier_type_id());
73
74 if (quan_.get_unknown_mask(elem))
75 {
76 for (std::size_t index_H=1; index_H < quan_.get_value_H_size() - 1; ++index_H)
77 {
78 const long unknown_index = quan_.get_unknown_index(elem, index_H);
79 if (unknown_index < 0)
80 continue;
81
82 const double energy_lower = std::max( (quan_.get_kinetic_energy(elem, index_H - 1) + quan_.get_kinetic_energy(elem, index_H)) / 2.0,
83 energy_start_);
84 const double energy_upper = std::min( (quan_.get_kinetic_energy(elem, index_H + 1) + quan_.get_kinetic_energy(elem, index_H)) / 2.0,
85 energy_end_);
86
87 if (energy_upper >= 0 && energy_lower < energy_upper)
88 {
89 double height = box_height(quan_, elem, index_H);
90 switch (conf_.she_discretization_type())
91 {
93 density += quan_.get_values(elem, index_H)[0] * averaged_density_of_states(quan_, dispersion, elem, index_H) * height;
94 break;
96 density += quan_.get_values(elem, index_H)[0] * height;
97 break;
98 default: throw std::runtime_error("carrier_density_wrapper_by_reference::operator(): Unknown SHE discretization type!");
99 }
100 }
101 } // for index_H
102 }
103 else // If there's a boundary condition
104 {
105 for (std::size_t index_H=1; index_H < quan_.get_value_H_size() - 1; ++index_H)
106 {
107 const double energy_lower = std::max( (quan_.get_kinetic_energy(elem, index_H - 1) + quan_.get_kinetic_energy(elem, index_H)) / 2.0,
108 energy_start_);
109 const double energy_upper = std::min( (quan_.get_kinetic_energy(elem, index_H + 1) + quan_.get_kinetic_energy(elem, index_H)) / 2.0,
110 energy_end_);
111 if (energy_upper >= 0 && energy_lower < energy_upper)
112 {
113 double height = box_height(quan_, elem, index_H);
114 switch (conf_.she_discretization_type())
115 {
117 density += quan_.get_boundary_value(elem, index_H) * averaged_density_of_states(quan_, dispersion, elem, index_H) * height;
118 break;
120 density += quan_.get_boundary_value(elem, index_H) * height;
121 break;
122 default: throw std::runtime_error("carrier_density_wrapper_by_reference::operator(): Unknown SHE discretization type!");
123 }
124 }
125 } // for index_H
126 }
127
128 if ( viennashe::util::is_Inf(density) ) viennashe::log::warn() << "* carrier_density_wrapper()::operator(): WARNING: density is inf " << elem << std::endl;
129 if ( viennashe::util::is_NaN(density) ) viennashe::log::warn() << "* carrier_density_wrapper()::operator(): WARNING: density is nan " << elem << std::endl;
130
131 return density * dispersion.symmetry_factor();
132 }
133
134 private:
135 viennashe::config conf_; // NO REFERENCE!
136 SHEQuantity const & quan_;
137 double energy_start_;
138 double energy_end_;
139 };
140 } // namespace detail
141
143 template <typename SHEQuantity>
145 {
148
149 public:
150
152
154 SHEQuantity const & quan,
155 double energy_start = 0.0,
156 double energy_end = 1.0)
157 : quan_(quan /* COPY ALL VALUES */),
158 density_impl_by_ref_(conf, quan_ /* REFERENCE COPIED VALUES */, energy_start, energy_end) { }
159
161 : quan_(o.quan_ /* COPY ALL VALUES */),
162 density_impl_by_ref_(o.density_impl_by_ref_, quan_ /* REFERENCE COPIED VALUES */)
163 { }
164
165 template <typename ElementType>
166 value_type operator()(ElementType const & elem) const
167 {
168 return this->density_impl_by_ref_(elem);
169 }
170
171 private:
172 SHEQuantity quan_;
173 carrier_density_by_ref_type density_impl_by_ref_;
174 };
175
176
178 template <typename DeviceType,
179 typename SHEQuantity,
180 typename VectorType>
181 void compute_carrier_density_vector(DeviceType const & device,
182 SHEQuantity const & quan,
184 VectorType & carrier)
185 {
186 typedef typename DeviceType::mesh_type MeshType;
187 typedef typename viennagrid::result_of::const_cell_range<MeshType>::type CellContainer;
188 typedef typename viennagrid::result_of::iterator<CellContainer>::type CellIterator;
189
190 MeshType const & mesh = device.mesh();
191
192 carrier_density_wrapper<SHEQuantity> carrier_dens(quan, dispersion);
193
194 //carrier concentration is computed for even unknowns, which are associated with vertices:
195 CellContainer cells(mesh);
196 for (CellIterator cit = cells.begin();
197 cit != cells.end();
198 ++cit)
199 {
200 long carrier_index = quan.get_unknown_index(*cit);
201
202 //only compute in the interior
203 if (carrier_index >= 0)
204 carrier[carrier_index] = carrier_dens(*cit);
205 }
206 }
207
208
218 template <typename DeviceType,
219 typename SHEQuantity,
220 typename ContainerType>
222 SHEQuantity const & quan,
224 ContainerType & container,
225 double energy_start = 0.0,
226 double energy_end = 1.0)
227 {
228 carrier_density_wrapper<SHEQuantity> wrapper(quan, dispersion, energy_start, energy_end);
229
231 }
232
233
234
235 } //namespace she
236} //namespace viennashe
237
238#endif
Defines a set of checker functors for micro-tests within ViennaSHE.
The main SHE configuration class. To be adjusted by the user for his/her needs.
Definition: config.hpp:124
viennashe::physics::dispersion_proxy dispersion_relation(viennashe::carrier_type_id ctype) const
Returns the dispersion relation for electrons.
Definition: config.hpp:266
she_discretization_type_id she_discretization_type() const
Definition: config.hpp:360
MeshT const & mesh() const
Returns the underlying mesh.
Definition: device.hpp:145
Defines the physical properties of a device, e.g. doping. This is the implementation for 2d and highe...
Definition: device.hpp:818
A proxy object for a dispersion relation. Does NOT take ownership of the provided pointer!
Definition: dispersion.hpp:69
An accessor for the carrier density in the device.
carrier_density_wrapper(carrier_density_wrapper const &o)
carrier_density_by_ref_type::value_type value_type
value_type operator()(ElementType const &elem) const
carrier_density_wrapper(viennashe::config const &conf, SHEQuantity const &quan, double energy_start=0.0, double energy_end=1.0)
An accessor for the carrier density in the device by reference
carrier_density_wrapper_by_reference(carrier_density_wrapper_by_reference const &o, SHEQuantity const &quan)
value_type operator()(ElementType const &elem) const
carrier_density_wrapper_by_reference(carrier_density_wrapper_by_reference const &o)
carrier_density_wrapper_by_reference(viennashe::config const &conf, SHEQuantity const &quan, double energy_start=0.0, double energy_end=1.0)
The SHE configuration class is defined here.
Contains forward declarations and definition of small classes that must be defined at an early stage.
Implementation of numerical integration routines.
Writer for arbitrary macroscopic quantities.
logger< true > warn()
Used to log warnings. The logging level is logWARNING.
Definition: log.hpp:303
double averaged_density_of_states(SHEQuantity const &quan, viennashe::config::dispersion_relation_type const &dispersion, CellFacetType const &cell_facet, std::size_t index_H)
Returns the density of states around a vertex or an edge at total energy specified by index_H....
void compute_carrier_density_vector(DeviceType const &device, SHEQuantity const &quan, viennashe::config::dispersion_relation_type const &dispersion, VectorType &carrier)
Computes the carrier density in the device and writes the result to a vector. Used during Gummel iter...
double box_height(SHEQuantity const &quan, ElementType const &elem, std::size_t index_H)
Returns the height of the control box with respect to energy at the provided element (vertex or edge)...
void write_carrier_density_to_container(DeviceType const &device, SHEQuantity const &quan, viennashe::config::dispersion_relation_type const &dispersion, ContainerType &container, double energy_start=0.0, double energy_end=1.0)
Convenience function for writing the average expansion order to the container provided.
bool is_Inf(const ValueType &val)
Checks if a value of type ValueType is Inf using (value - value) != 0.
Definition: checks.hpp:101
bool is_NaN(const ValueType &val)
Checks if a value of type ValueType is NaN using value != value.
Definition: checks.hpp:97
The main ViennaSHE namespace. All functionality resides inside this namespace.
Definition: accessors.hpp:40
@ SHE_DISCRETIZATION_EVEN_ODD_ORDER_GENERALIZED_DF
Definition: forwards.h:146
@ SHE_DISCRETIZATION_EVEN_ODD_ORDER_DF
Definition: forwards.h:145
void write_macroscopic_quantity_to_container(DeviceType const &device, MacroscopicQuantityAccessor const &quantity, ContainerType &quantity_container)
Writes the provided macroscopic quantity to the container provided.
Definition: macroscopic.hpp:46
Provides a number of fundamental constants. All constants in SI units.
Implementation of spherical harmonics plus helper functions.