PLaSK library
|
Solver performing calculations in 2D Cartesian or Cylindrical space using finite element method. More...
#include <solvers/electrical/ddm2d/ddm2d.hpp>
Public Member Functions | |
double | compute (unsigned loops=1) |
Run drift_diffusion calculations. | |
double | findEnergyLevels () |
Find energy levels - TEST. | |
bool | checkWell (std::string _carrier) |
Check if the carriers are confined. | |
double | integrateCurrent (size_t vindex, bool onlyactive=false) |
Integrate vertical total current at certain level. | |
double | getTotalCurrent (size_t nact=0) |
Integrate vertical total current flowing vertically through active region. | |
void | loadConfiguration (XMLReader &source, Manager &manager) override |
Compute total electrostatic energy stored in the structure. | |
DriftDiffusionModel2DSolver (const std::string &name="") | |
std::string | getClassName () const override |
Get name of solver. | |
~DriftDiffusionModel2DSolver () | |
double | integrateCurrent (size_t vindex, bool onlyactive) |
double | integrateCurrent (size_t vindex, bool onlyactive) |
std::string | getClassName () const |
Get name of solver. | |
std::string | getClassName () const |
Get name of solver. | |
Public Member Functions inherited from plask::FemSolverWithMesh< Geometry2DType, RectangularMesh< 2 > > | |
FemSolverWithMesh (const std::string &name="") | |
bool | parseFemConfiguration (XMLReader &reader, Manager &manager) |
FemMatrix * | getMatrix () |
FemMatrix * | getMatrix () |
Public Member Functions inherited from plask::SolverWithMesh< SpaceT, MeshT > | |
SolverWithMesh (const std::string &name="") | |
~SolverWithMesh () | |
void | loadConfiguration (XMLReader &source, Manager &manager) override |
Load configuration from given source . | |
void | parseStandardConfiguration (XMLReader &source, Manager &manager, const std::string &expected_msg="solver configuration element") |
virtual void | onMeshChange (const typename MeshT::Event &PLASK_UNUSED(evt)) |
This method is called just after the mesh has been changed. | |
void | onGeometryChange (const Geometry::Event &PLASK_UNUSED(evt)) override |
This method is called when the geometry is changed. | |
MeshT & | meshRef () const |
Get current module mesh. | |
shared_ptr< MeshT > | getMesh () const |
Get current solver mesh. | |
void | setMesh (const shared_ptr< MeshT > &mesh) |
Set new mesh for the solver. | |
void | setMesh (shared_ptr< MeshGeneratorD< MeshT::DIM > > generator) |
Set new mesh got from generator. | |
Public Member Functions inherited from plask::SolverOver< SpaceT > | |
SolverOver (const std::string &name="") | |
~SolverOver () | |
void | parseStandardConfiguration (XMLReader &source, Manager &manager, const std::string &expected_msg="solver configuration element") |
shared_ptr< SpaceT > | getGeometry () const |
Get current solver geometry space. | |
void | setGeometry (const shared_ptr< SpaceT > &geometry) |
Set new geometry for the solver. | |
Public Member Functions inherited from plask::Solver | |
bool | initCalculation () |
This should be called on beginning of each calculation method to ensure that solver will be initialized. | |
Solver (const std::string &name="") | |
Construct uninitialized solver. | |
virtual | ~Solver () |
Virtual destructor (for subclassing). Do nothing. | |
void | parseStandardConfiguration (XMLReader &source, Manager &manager, const std::string &expected_msg="solver configuration element") |
Load standard configuration (geometry, mesh) tags from source . | |
bool | isInitialized () |
Check if solver is already initialized. | |
void | invalidate () |
This method should be and is called if something important was changed: calculation space, mesh, etc. | |
std::string | getId () const |
Get solver id. | |
std::string | getName () const |
virtual std::string | getClassDescription () const |
Get a description of this solver. | |
template<typename ArgT = double, typename ValT = double> | |
DataLog< ArgT, ValT > | dataLog (const std::string &chart_name, const std::string &axis_arg_name, const std::string &axis_val_name) |
template<typename ArgT = double, typename ValT = double> | |
DataLog< ArgT, ValT > | dataLog (const std::string &axis_arg_name, const std::string &axis_val_name) |
template<typename ... Args> | |
void | writelog (LogLevel level, std::string msg, Args &&... params) const |
Log a message for this solver. | |
Public Attributes | |
double | maxerr |
bands have to be up-shifted - we want only positive values of energy levels; unit: eV | |
BoundaryConditions< RectangularMesh< 2 >::Boundary, double > | voltage_boundary |
Boundary condition. | |
ProviderFor< Potential, Geometry2DType >::Delegate | outPotential |
ProviderFor< FermiLevels, Geometry2DType >::Delegate | outFermiLevels |
ProviderFor< BandEdges, Geometry2DType >::Delegate | outBandEdges |
ProviderFor< CurrentDensity, Geometry2DType >::Delegate | outCurrentDensityForElectrons |
ProviderFor< CurrentDensity, Geometry2DType >::Delegate | outCurrentDensityForHoles |
ProviderFor< CarriersConcentration, Geometry2DType >::Delegate | outCarriersConcentration |
ProviderFor< Heat, Geometry2DType >::Delegate | outHeat |
ReceiverFor< Temperature, Geometry2DType > | inTemperature |
bool | mRsrh |
SRH recombination is taken into account. | |
bool | mRrad |
radiative recombination is taken into account | |
bool | mRaug |
Auger recombination is taken into account. | |
bool | mPol |
polarization (GaN is the substrate) | |
bool | mFullIon |
dopant ionization = 100% | |
double | mSchottkyP |
Schottky barrier for p-type contact (eV) | |
double | mSchottkyN |
Schottky barrier for n-type contact (eV) | |
double | maxerrPsiI |
Maximum estimated error for initial potential during all iterations (useful for single calculations managed by external python script) | |
double | maxerrPsi0 |
Maximum estimated error for potential at U = 0 V during all iterations (useful for single calculations managed by external python script) | |
double | maxerrPsi |
Maximum estimated error for potential during all iterations (useful for single calculations managed by external python script) | |
double | maxerrFn |
Maximum estimated error for quasi-Fermi energy level for electrons during all iterations (useful for single calculations managed by external python script) | |
double | maxerrFp |
Maximum estimated error for quasi-Fermi energy level for holes during all iterations (useful for single calculations managed by external python script) | |
size_t | loopsPsiI |
Loops limit for initial potential. | |
size_t | loopsPsi0 |
Loops limit for potential at U = 0 V. | |
size_t | loopsPsi |
Loops limit for potential. | |
size_t | loopsFn |
Loops limit for quasi-Fermi energy level for electrons. | |
size_t | loopsFp |
Loops limit for quasi-Fermi energy level for holes. | |
Public Attributes inherited from plask::FemSolverWithMesh< Geometry2DType, RectangularMesh< 2 > > | |
FemMatrixAlgorithm | algorithm |
Factorization algorithm to use. | |
IterativeMatrixParams | iter_params |
Parameters of iterative solver. | |
Protected Member Functions | |
void | onInitialize () override |
Initialize the solver. | |
void | onInvalidate () override |
Invalidate the data. | |
size_t | getActiveRegionMeshIndex (size_t actnum) const |
Get info on active region. | |
void | onMeshChange (const typename RectangularMesh< 2 >::Event &evt) override |
void | onGeometryChange (const Geometry::Event &evt) override |
This method is called when the geometry is changed. | |
void | computePsiI () |
Calculate initial potential for all elements. | |
size_t | isActive (const Vec< 2 > &point) const |
Return true if the specified point is at junction. | |
size_t | isActive (const RectangularMesh< 2 >::Element &element) const |
Return true if the specified element is a junction. | |
const LazyData< double > | getPotentials (shared_ptr< const MeshD< 2 > > dest_mesh, InterpolationMethod method) const |
const LazyData< double > | getFermiLevels (FermiLevels::EnumType what, shared_ptr< const MeshD< 2 > > dest_mesh, InterpolationMethod method) const |
const LazyData< double > | getBandEdges (BandEdges::EnumType what, shared_ptr< const MeshD< 2 > > dest_mesh, InterpolationMethod method) |
const LazyData< double > | getHeatDensities (shared_ptr< const MeshD< 2 > > dest_mesh, InterpolationMethod method) |
const LazyData< Vec< 2 > > | getCurrentDensitiesForElectrons (shared_ptr< const MeshD< 2 > > dest_mesh, InterpolationMethod method) |
const LazyData< Vec< 2 > > | getCurrentDensitiesForHoles (shared_ptr< const MeshD< 2 > > dest_mesh, InterpolationMethod method) |
const LazyData< double > | getCarriersConcentration (CarriersConcentration::EnumType what, shared_ptr< const MeshD< 2 > > dest_mesh, InterpolationMethod method) |
Protected Member Functions inherited from plask::SolverOver< SpaceT > | |
template<typename Boundary , typename ConditionT > | |
void | readBoundaryConditions (Manager &manager, XMLReader &reader, BoundaryConditions< Boundary, ConditionT > &dest) |
Read boundary conditions using information about the geometry of this solver. | |
Protected Attributes | |
std::size_t | size |
Number of columns in the main matrix. | |
double | mTx |
ambient temperature (K) | |
double | mEx |
energy (eV) | |
double | mNx |
maximal doping concentration (1/cm^3) | |
double | mEpsRx |
maximal dielectric constant (-) | |
double | mXx |
sometimes denoted as LD (um) | |
double | mMix |
maximal mobility (cm^2/Vs) | |
double | mRx |
recombination parameter (1/(cm^3*s)) | |
double | mJx |
current density parameter (kA/cm2) | |
double | mAx |
radiative recombination coefficient (1/s) | |
double | mBx |
radiative recombination coefficient (cm^3/s) | |
double | mCx |
Auger recombination coefficient (cm^6/s) | |
double | mPx |
polarization (C/m^2) | |
double | dU |
default voltage step (V) | |
double | maxDelPsi0 |
maximal correction for initial potential calculations (V) | |
double | maxDelPsi |
maximal correction for potential calculations (V) | |
double | maxDelFn |
maximal correction for quasi-Fermi levels for electrons calculations (eV) | |
double | maxDelFp |
maximal correction for quasi-Fermi levels for holes calculations (eV) | |
Stat | stat |
carriers statistics | |
ContType | conttype |
type of contacts (ohmic/Schottky) | |
DataVector< double > | dveN |
Cached electron concentrations (size: elements) | |
DataVector< double > | dveP |
Cached hole concentrations (size: elements) | |
DataVector< double > | dvePsi |
Computed potentials (size: elements) | |
DataVector< double > | dveFnEta |
Computed exponents of quasi-Fermi levels for electrons (size: elements) | |
DataVector< double > | dveFpKsi |
Computed exponents of quasi-Fermi levels for holes (size: elements) | |
DataVector< double > | dvnPsi0 |
Computed potential for U=0V (size: nodes) | |
DataVector< double > | dvnPsi |
Computed potentials (size: nodes) | |
DataVector< double > | dvnFnEta |
Computed exponents of quasi-Fermi levels for electrons (size: nodes) | |
DataVector< double > | dvnFpKsi |
Computed exponents of quasi-Fermi levels for holes (size: nodes) | |
DataVector< Vec< 2, double > > | currentsN |
Computed current densities for electrons. | |
DataVector< Vec< 2, double > > | currentsP |
Computed current densities for holes. | |
DataVector< double > | heats |
Computed and cached heat source densities. | |
bool | needPsi0 |
Flag indicating if we need to compute initial potential;. | |
double | T0 |
Temperature used for compiting level estimates. | |
bool | strained |
Consider strain in QW? | |
Protected Attributes inherited from plask::SolverWithMesh< SpaceT, MeshT > | |
shared_ptr< MeshT > | mesh |
Mesh over which the calculations are performed. | |
boost::signals2::connection | mesh_signal_connection |
Connection of mesh to onMeshChange method, see http://www.boost.org/doc/libs/1_55_0/doc/html/signals2/tutorial.html#idp204830936. | |
Protected Attributes inherited from plask::SolverOver< SpaceT > | |
shared_ptr< SpaceT > | geometry |
Space in which the calculations are performed. | |
Protected Attributes inherited from plask::Solver | |
bool | initialized |
true only if solver is initialized | |
Additional Inherited Members | |
Public Types inherited from plask::SolverWithMesh< SpaceT, MeshT > | |
typedef MeshT | MeshType |
Type of the mesh for this solver. | |
Public Types inherited from plask::SolverOver< SpaceT > | |
typedef SpaceT | SpaceType |
of the space for this solver | |
Solver performing calculations in 2D Cartesian or Cylindrical space using finite element method.
plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::DriftDiffusionModel2DSolver | ( | const std::string & | name = "" | ) |
bool plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::checkWell | ( | std::string | _carrier | ) |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::compute | ( | unsigned | loops = 1 | ) |
Run drift_diffusion calculations.
|
protected |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::findEnergyLevels | ( | ) |
Find energy levels - TEST.
hb*hb/(2m), unit: eV*nm*nm, 10^9 is introduced to change m into nm
bands have to be up-shifted - we want only positive values of energy levels; unit: (eV)
TODO
true when electrons (el) are confined
true when heavy holes (hh) are confined
true when light holes (lh) are confined
1/(dz*dz), unit: (1/nm^2)
1/(dz), unit: (1/nm)
the order of the small matrix for central-node for CB
the order of the matrix for CB
N - Hc size
setting Hc = zeros
nn - number of nodes
centre of the left element
centre of the right element
left/central/right 1x1 local matrix
11
Re - real part., Im - imag part.
no energy levels for electrons
lev_el - vector with energy levels for electrons
n_lev_el - number of energy levels of electrons
we have to find out if this level corresponds to el
"1" in both cases because only el are considered here
sorting electron energy levels
|
protected |
|
protected |
|
protected |
|
virtual |
|
virtual |
|
overridevirtual |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::getTotalCurrent | ( | size_t | nact = 0 | ) |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DCartesian >::integrateCurrent | ( | size_t | vindex, |
bool | onlyactive | ||
) |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DCylindrical >::integrateCurrent | ( | size_t | vindex, |
bool | onlyactive | ||
) |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::integrateCurrent | ( | size_t | vindex, |
bool | onlyactive = false |
||
) |
Integrate vertical total current at certain level.
vindex | vertical index of the element mesh to perform integration at |
onlyactive | if true only current in the active region is considered |
|
inlineprotected |
|
inlineprotected |
|
overridevirtual |
Compute total electrostatic energy stored in the structure.
Reimplemented from plask::Solver.
|
inlineoverrideprotectedvirtual |
This method is called when the geometry is changed.
It just calls invalidate(); but subclasses can customize it.
evt | information about geometry changes |
Reimplemented from plask::SolverOver< SpaceT >.
|
overrideprotectedvirtual |
|
overrideprotectedvirtual |
|
inlineoverrideprotected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
ReceiverFor<Temperature, Geometry2DType> plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::inTemperature |
size_t plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::loopsFn |
size_t plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::loopsFp |
size_t plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::loopsPsi |
size_t plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::loopsPsi0 |
size_t plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::loopsPsiI |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::maxerrFn |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::maxerrFp |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::maxerrPsi |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::maxerrPsi0 |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::maxerrPsiI |
|
protected |
|
protected |
|
protected |
|
protected |
bool plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::mFullIon |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::mSchottkyN |
double plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::mSchottkyP |
|
protected |
|
protected |
|
protected |
ProviderFor<BandEdges,Geometry2DType>::Delegate plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::outBandEdges |
ProviderFor<CarriersConcentration,Geometry2DType>::Delegate plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::outCarriersConcentration |
ProviderFor<CurrentDensity,Geometry2DType>::Delegate plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::outCurrentDensityForElectrons |
ProviderFor<CurrentDensity,Geometry2DType>::Delegate plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::outCurrentDensityForHoles |
ProviderFor<FermiLevels,Geometry2DType>::Delegate plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::outFermiLevels |
ProviderFor<Heat,Geometry2DType>::Delegate plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::outHeat |
ProviderFor<Potential,Geometry2DType>::Delegate plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::outPotential |
|
protected |
|
protected |
|
protected |
|
protected |
Temperature used for compiting level estimates.
BoundaryConditions<RectangularMesh<2>::Boundary, double> plask::electrical::drift_diffusion::DriftDiffusionModel2DSolver< Geometry2DType >::voltage_boundary |