|
| |
Linear current equation
In this chapter it is described how to set up the current matrix for the
linear current equation:
The following topics are addressed:
- How current-clusters are defined and set up.
- How the current matrix is assembled.
- How boundary conditiond are implemented.
- How the current matrix is stored.
- Special treatment of current in 1D.
- How to calculate the matrix-vector product.
Current clusters
Before the current matrix can be created it must
be defined
- where the current equation is solved.
- which boundary conditions on domain/inner
boundaries to take.
The relevant information is written to
MODULE current:
[first half]
Current clusters are - as Poisson contact clusters - the Poisson contact regions defined via
the input system. Poisson contact regions are a set of points on the material
grid. All adjoining points
on the physical grid are assigned the corresponding boundary condition. The
Poisson contact region is NOT the area where the Poisson equation has to be
solved, the Poisson contact region defines the contacts specified via the keyword
$poisson-boundary-conditions.
MODULE
current
num_current_clusters ... Number of
current clusters
current_clusterM(1..num_current_clusters) has TYPE
TYPE :: current_cluster
INTEGER :: num_points Number of points in that cluster
INTEGER,DIMENSION(1:num_points),POINTER :: num_guenther_points
Number of points on
material grid which belong to that cluster.
INTEGER :: kind_point = 1 -->
Dirichlet condition for Fermi level
Fermi level set to value [V]
= 2 -->
set to Dirichlet value,
Fermi equilibrium value in [V]
= 3 -->
set to Neumann value, [A/m2] for current
REAL(8),DIMENSION(3) :: valueV ... see
kind_point
Note: The field is a vector E =
(Ex,Ey,Ez) = Ex * x + Ey *
y +
Ez * z.
MODULE
current
The Poisson assembler takes all 8 (3D) octants for each
point.
Wherever there is a transition from an octant on a
Poisson cluster to an octant outside a Poisson cluster and a field is specified
for that cluster, the corresponding component of the field is taken.
If any of the 8 octants has the option kind_point=2,3,4
(i.e. any kind
of Dirichlet) the assembler ignores all fields and
sets the potential on that
point to the specified Dirichlet value.
For kind_point=4 --> the band edge of
the adjoining octant is taken.
Note: This option is used for
Schottky contacts, i.e. at a contact/semiconductor
interface. So one has to take the band edge in the semiconductor, not in
the contact.
If
colliding potentials are specified the cluster with the lowest value
is valid.
The same procedure for boundary points on boundary
of simulation domain:
currentbnd_x3DM(1..2,j,k,4): currentbnd_x3DM(1,j,k)
--> i=1
currentbnd_x3DM(2,j,k) --> i=n_1
1..4 --> interfaces 9,10,11,12
currentbnd_y3DM(1..2,i,k,4): currentbnd_x3DM(1,i,k)
--> j=1
currentbnd_x3DM(2,i,k) --> j=n_2
1..4 --> interfaces 5,6,7,8
currentbnd_z3DM(1..2,i,j,4):
currentbnd_x3DM(1,i,j) --> k=1
currentbnd_x3DM(2,i,j) --> k=n_1
1..4 --> interfaces 1,2,3,4
Note: Numbering of interfaces (1..12)
interface 1 <--> between octant 1 and 5
interface 2 <--> between octant 2 and 6
interface 3 <--> between octant 3 and 7
interface 4 <--> between octant 4 and 8
interface 5 <--> between octant 1 and 3
interface 6 <--> between octant 2 and 4
interface 7 <--> between octant 5 and 7
interface 8 <--> between octant 6 and 8
interface 9 <--> between octant 1 and 2
interface 10 <--> between octant 3 and 4
interface 11 <--> between octant 5 and 6
interface 12 <--> between octant 7 and 8
In 2D:
Note: Numbering of interfaces (1..4)
interface 1 <--> between
octant 1 and 2
interface 2 <--> between
octant 1 and 3
3
interface 3 <--> between
octant 3 and 4
2 + 4
interface 4 <--> between
octant 2 and 4 1
= 1 -->
Dirichlet condition for Fermi level, Fermi level set to value [V]
= 2 -->
set to Dirichlet value, Fermi
equilibrium value in [V]
= 3 -->
set to Neumann value, [A/m2]
for current
bndval_x3DM(1..2,j,k,4)
bndval_y3DM(1..2,i,k,4)
bndval_z3DM(1..2,i,j,4)
... field, if
currentbnd_.3DM=1
(directed orthonormal to the surface)
potential currentbnd_.3DM=2 ("."
stands for x, y or z)
Note: Up to now, on simulation domain boundary only field zero is implemented.
Must be allocated and defined before call of
currentmatrix3D
and is deallocated in
deallocate_current.
The current clusters and contact models are
provided in the modules:
mod_num_of_ohmics
mod_num_of_schottkys
mod_num_of_dirichlets
mod_num_of_neumanns
mod_num_of_mixed
mod_apply_dirichlet_to_cluster
mod_dirichlet_grid_points
mod_potential_for_dirichlet
mod_apply_ohmic_to_cluster
mod_ohmic_grid_points
mod_voltage_for_ohmic
mod_apply_schottky_to_cluster
mod_schottky_grid_points
mod_voltage_for_schottky
mod_barrier_for_schottky
See there for more details.
With that information subroutine
define_current_clusters defines the current problem in MODULE
current:
ohmic, voltage
controlled --> |
When voltage is applied, Dirichlet condition
(for potential and Fermi levels) with
pot = pot_built_in + V_applied
fermi = fermi_equilibrium - V_applied. |
ohmic, current
controlled --> |
not yet
implemented |
Schottky,
voltage controlled --> |
When voltage is applied, Dirichlet condition
(for potential and Fermi levels) with
pot = pot_built_in + V_applied
fermi = fermi_equilibrium - V_applied. |
Schottky,
current controlled --> |
not yet
implemented |
Dirichlet
--> |
not yet
implemented |
Neumann
--> |
not yet
implemented |
mixed, i.e.
relation between Dirichlet
and Neumann (Cauchy)
--> |
not yet
implemented |
Note: The grid points inside a current
cluster (contact) are exempted from simulation (i.e. set to Dirichlet which is
of no importance).
Setup of current matrix
After the Poisson clusters have been defined by
subroutine
define_current_clusters the current matrix is set up in
SUBROUTINE currentmatrix3D
(currentmatrix2D,currentmatrix1D).
The current solver allows one to exempt parts of the
simulation domain from simulation (contacts). On these internal surfaces Neumann boundary
conditions can be specified (zero by
default).
The corresponding information is provided by
SUBROUTINE
irregular_boundary('Current',i,j,k,irr_bcC,which_octantsVL,fieldM):
Note: Up to now only Dirichlet boundary conditions on contacts are used.
SUBROUTINE
irregular_boundary(i,j,k,irr_bcC,which_octantsVL,fieldM)
(i,j,k) point on physical grid
[Input]
which_octantsVL(1..8) .TRUE. --->
octant
within simulation region [Output]
.FALSE. ---> octant outside simulation region, i.e.
inside contact
irr_bcC = 'in' If all octants within
simulation region.
'ot' If all octants outside simulation region.
'on' If
at least one octant within simulation region, but not all.
curM(1..8,1..3) : curM(n,1) field component
in x-direction for box n
curM(n,2) field component in y-direction for box n
curM(n,3)
field component in z-direction for box n
[units in
A/m2]
NOTE: Up to now only Dirichlet conditions on contacts are used. So curM
is zero.
Finally, the program in first SETUP always sets
Neumann boundary conditions which are set to zero if none are defined. Later on, Dirichlet boundary
conditions are set (i.e. Neumann are overwritten by Dirichlet).
The fields on the domain boundaries can be obtained
via
boundary_value:
FUNCTION
boundary_value(ch,i,j,k,whichinterface) RESULT(bv)
This may just be called for (i,j,k) on the domain boundary, otherwise stop.
Input:
(i,j,k) on physical grid
ch = 'x' or 'y' or 'z' (direction of differential)
called at any setup of current matrix
whichinterface=1..4 --> interfaces
9,10,11,12 for 'x'
5, 6, 7, 8
for 'y'
1, 2, 3, 4
for 'z'
For whichinterface the same convention as for
currentbnd_x3DM... is
.TRUE.
Output:
If Neumann boundary condidtions are valid, the current is provided.
Up to now: not used --> set to zero
Otherwise: set to zero
Note: In the first setup always Neumann conditions are implemented.
They are overwritten with Dirichlet values (if specified) later on.
The mobility is provided by subroutine
get_mobility3D(get_mobility2D).
The mobility is needed on a point on the material grid. The specification is
done as described above.
SUBROUTINE
get_mobility3D(i1,j1,k1,oct1,i2,j2,k2,oct2,chargeC,mobility)
Input:
(i1,j1,k1) point on physical
grid, oct1 octant(1..8)
(i2,j2,k2) point on physical
grid, oct2 octant(1..8)
3 |
4
3 | 4
|
|
---(i1,j1)--- ---(i2,j2)---
|
|
1 |
2
1 | 2
chargeC = 'el' or 'hl'
Output:
Mobility between the specified two octants in [cm2/Vs]
Note: They must be adjoining points and octant numbering must be consistent.
The recombination is provided by subroutine
get_net_recombination.
It refers to a point of the material grid.
SUBROUTINE
get_net_recombination(i,j,k,oct,chargeC,net_recomb)
Input:
(i,j,k) point on Physical grid
oct ... 1..8 number of octant
3 | 4
|
----(i,j)----
|
1 | 2
chargeC = 'el' or 'hl'
Output:
net_recomb ... net recombination rate in [e/m3/s]
Comprises recombination and generation terms.
- is positive for effective generation
- is negative for effective recombination
Note:
This routine is used for option 'DriftDiffAll'.
This means that one quasi-Fermi level is assumed for all carriers (which is the
only option up to now): One Fermi level for all electrons and another Fermi
level for all holes.
The computational molecule is assembled by first calculating the
computational molecules for all eight octants and adding them later on.
Subroutine set_cur_differential3D
(set_cur_differential2D,
set_cur_differential1D,
FUNCTION differential
(1D only))
provides the computational molecule for a specified octant.
SUBROUTINE set_cur_differential3D(i,j,k,oct,chargeC,modeC,starV)
Defines differential for current equation for octant oct of point
(i,j,k).
Input:
i,j,k point on physical grid
oct 1...8 number of octant
modeC
NOTE: Numbering of octants
3 | 4
|
-----(i,j,k-)-----
|
1 | 2
7 | 8
|
-----(i,j,k+)-----
|
5 | 6
chargeC = 'el' or 'hl'
Output:
The computational molecule is provided in starV.
At first, each octant has its own numbering of points as specified below. For each octant there are four neighboring points that are needed for finite differences. Their numbering is as follows:
octant1 (i ,j ,k-1) = 1 (i-,j-,k-)
(i ,j-1,k ) = 2
(i-1,j ,k ) = 3
(i ,j ,k ) = 4
octant2 (i ,j ,k-1) = 5 (i+,j-,k-)
(i ,j-1,k ) = 6
(i+1,j ,k ) = 7
(i ,j ,k ) = 8
octant3 (i ,j ,k-1) = 9 (i-,j+,k-)
(i ,j+1,k ) =
10
(i-1,j ,k ) = 11
(i ,j ,k ) = 12
octant4 (i ,j ,k-1) = 13 (i+,j+,k-)
(i ,j+1,k ) =
14
(i+1,j ,k ) =
15
(i ,j ,k ) = 16
octant5 (i ,j ,k+1) = 17 (i-,j-,k+)
(i ,j-1,k ) = 18
(i-1,j , k) = 19
(i ,j ,k ) = 20
octant6 (i ,j ,k+1) = 21 (i+,j-,k+)
(i ,j-1,k ) = 22
(i+1,j ,k ) = 23
(i ,j ,k ) = 24
octant7 (i ,j ,k+1) = 25 (i-,j+,k+)
(i ,j+1 ,k) = 26
(i-1,j ,k ) = 27
(i ,j ,k ) = 28
octant8 (i ,j ,k+1) = 29 (i+,j+,k+)
(i ,j+1,k ) = 30
(i+1,j ,k ) = 31
(i ,j ,k ) = 32
right handside: octant1 = 33
octant2 = 34
octant3 = 35
octant4 = 36
octant5 = 37
octant6 = 38
octant7 = 39
octant8 = 40
starV(1:40)
if oct = 1 --> starV( 1: 4) and starV(33)
is defined, rest=0
= 2 --> starV( 5: 8) and starV(34)
is
defined, rest=0
= 3 --> starV( 9:12) and
starV(35) is defined, rest=0
= 4 --> starV(13:16) and
starV(36) is defined, rest=0
= 5 --> starV(17:20) and starV(37)
is
defined, rest=0
= 6 --> starV(21:24) and starV(38)
is
defined, rest=0
= 7 --> starV(25:28) and
starV(39) is defined, rest=0
= 8 --> starV(29:32) and
starV(40) is defined, rest=0
The subroutine is called for each octant and the corresponding eight parts of
the computational molecule are defined. Later on they are added to build the final computational molecule.
Discretizes the differential operator
- int [div j] = - surface_integral j dA ==
rhs
where j = mobility * dens * grad(EF)
[A/m2]
rhs = int recombination
Note: In this first version we simply use the density for dens.
dens > 0, C/m3
Units:
mobility -->
m2/Vs
dens -->
C/m3 = As/m3
grad E_F -->
V/m --> j = A/m2
Note: EF is given in [V].
The matrix which is created at the setup is
stored in MODULE current:
diagV ..
DIMENSION(1:lengvecV)
contains diagonal elements of current matrix
matrixindicesM ..
DIMENSION(1:nummaelemente,1:2)
contains for nonvanishing offdiagonal matrix elements row and column
index
Because matrix is symmetric, only upper (triangle) half is stored.
matrixwerteV .. DIMENSION(1:nummaelemente)
contains corresponding values
matrixposV(i) ..
DIMENSION(1..lengvecV)
matrix_pos(n) points to the position, where the nonvanishing elements of row n
(upper triangle) are placed in the arrays.
matrixindicesM and matrixwerteV
If the line doesn't have any offdiagonal elements it points to the next entry
unequal to zero.
nummaelemente .. number of
offdiagonal matrix elements
bV ..
DIMENSION(1:lengvecV)
contains the complete right hand side of the current matrix equation
The setup is done in SUBROUTINE
currentmatrix3D (currentmatrix2D,currentmatrix1D):
The program loops through all points:
1) CALL WhichKindOfBoundaryPoint3D(ch,i,j,k)
[internal subroutine, in 2D/1D this distinction is done explicitly]
-->
ch ='in', if (i,j,k) is within domain.
otherwise --> ch provides information about the position of the boundary
point (see source code).
The
variables beteiligtL_V and randL_V provide
information for box integration:
beteiligtL_V(i)=.TRUE. <--> octant i is within simulation region, so box
integration is performed over it.
randL_V(1..3): says if flux through surface in x-, y-, z-direction has to
be considered in box integration.
2) CALL
irregular_boundary provides for (i,j,k) in irr_bcC the following information:
irr_bcC
= 'ot' <--> (i,j,k) lies in a
region which is exempted from simulation (i.e. a Poisson
contact region)
--> set all beteiligtL_V=0
= 'in' <--> (i,j,k) lies
within simulation region (i.e. not in a Poisson contact region)
--> do not alter beteiligtL_V
= 'on' <--> (i,j,k) lies on
internal interface separating
simulation / nonsimulation region (Poisson contact region)
--> set the corresponding octants of beteiligtL_V to 0
3) Now the setup of matrix elements for that
point:
a) if
outer point (irr_bcC=='ot') --> CALL outer_point,
writes for (i,j,k) a Dirichlet condition (=0)
b) if
inner point which is no multiple point and does not adjoin a multiple point and
not to a Poisson contact region:
the differential is defined explicitly: CALL standardfall3D
c) else: perform box integration over all boxes specified by
beteiligtL_V(:)=.TRUE.:
CALL CalculateCoefficient3D
4) Now allocate and define
ska_bnd_x3DM, ska_bnd_y3DM, ska_bnd_z3DM.
Technical details of computational molecule
First the computational molecule is extracted for each of the eight octants.
See description of subroutine set_cur_differential3D
for further information.
Finally the computational molecules are added to get the final computational
molecule:
1
2
3 4
5 6
7 |
8
---------------------------------------------------------------------------------------------------------------------
(i,j,k-1)
(i,j-1,k) (i-1,j,k) (i,j,k)
(i+1,j,k) (i,j+1,k) (i,j,k+1) |
right hand side
The column l=1,..,8 corresponds to koeffsum_V(l):
!--------------
! Addition of boxes
!--------------
koeffsum_V(1) = koeff_M( 1) + koeff_M( 5) + koeff_M( 9) + koeff_M(13)
koeffsum_V(2) = koeff_M( 2) + koeff_M( 6) + koeff_M(18) + koeff_M(22)
koeffsum_V(3) = koeff_M( 3) + koeff_M(11) + koeff_M(19) + koeff_M(27)
koeffsum_V(4) = koeff_M( 4) + koeff_M( 8) + koeff_M(12) + koeff_M(16) + &
koeff_M(20) + koeff_M(24) + koeff_M(28) + koeff_M(32)
koeffsum_V(5) = koeff_M( 7) + koeff_M(15) + koeff_M(23) + koeff_M(31)
koeffsum_V(6) = koeff_M(10) + koeff_M(14) + koeff_M(26) + koeff_M(30)
koeffsum_V(7) = koeff_M(17) + koeff_M(21) + koeff_M(25) + koeff_M(29)
koeffsum_V(8) = koeff_M(33) + koeff_M(34) + koeff_M(35) + koeff_M(36) + &
koeff_M(37) + koeff_M(38) + koeff_M(39) + koeff_M(40)
Now the area of internal interfaces are defined (only for multiple
points):
The matrix elements are written (just for
the upper part of the matrix (symmetric) in CreateMatrix3D:
Creates the elements for single points and the
first line of multiple points.
Note: For multiple points the balance
from integration is fully contained in the first line.
In
Poisson equation the potential is assumed to be continuous.
The vector of the Poisson equation is defined on Stefan's grid (1..lengvecV)
which contains
multiple point; so additional lines are added to account for that.
Since the matrix must be symmetric and positive definite, we use following
scheme:
Generate multiple points
Double points:
The matrix is supplemented by an identity relation for the second point. The
following structure is obtained:
2nd line: -a
+ a = 0
Matrix must be symmetric, thus add
+ a
and -a
to 1st line.
-1 0 ... 0 -1
4+a
-a -1 0
0 ... 0 -1 ...
= Q
... 0
0 -a
a 0 0 0 ... ...
= 0
-1 ...
0 0 -1 0
4 -1 0 ... 0 -1
... = Q
Eightfold point:
-1 ... -1 0 ... 0 -1
6+a -a
0 0 0 0 0 0 -1
... -1 ... -1 ... = Q
-a 7a
-a -a -a -a -a -a
= 0
0 -a a
0 0 0 0 0
= 0
0 -a 0
a 0 0 0
0
= 0
0 -a 0 0
a 0 0 0
= 0
0 -a 0 0
0 a 0 0
= 0
0 -a 0 0
0 0 a 0
= 0
0 -a 0 0
0 0 0 a
= 0
-1 ... -1 0 ... 0 -1
0 0 0 0 0 0
0 6 -1 .. -1 .. -1 .. = Q
(An alternative would be, but it didn't converge with conjugate gradient,
maybe later check again...)
-1 ... -1 0 ... 0 -1 6+7a -a -a -a -a -a -a -a -1 ... -1 ... -1
-a a
-a 0 a 0 0 0 0 0
-a 0 0 a 0 0 0 0
-a 0 0 0 a 0 0 0
-a 0 0 0 0 a 0 0
-a 0 0 0 0 0 a 0
-a 0 0 0 0 0 0 a
-1 ... -1 0 ... 0 -1 0 0 0 0 0 0 0 6 -1 .. -1 .. -1
Boundary Conditions (Dirichlet points):
After the matrix has been completely set up,
Neumann boundary conditions are overwritten with Dirichlet values:
There is an auxiliary subroutine DIRICHLET3D(i,j,k,dirval)
(2D/1D)
which
defines Dirichlet conditions for a certain point.
In detail: Matrix elements in that line are
overwritten with zero, diagonal element is set to a default value
and the
right hand side bV is defined
accordingly.
Note: Matrix must remain symmetric.
So if you overwrite one line to a Dirichlet condition, you also have
to change the line which shares the matrix element in order to restore symmetry.
SUBROUTINE
set_dirichlet_current
sets the
boundary conditions as specified in MODULE
current.
SUBROUTINE set_dirichlet_current
Sets Dirichlet points as specified in MODULE
current.
Current matrix must be defined.
Multiple points are all set to the same value.
The points with higher cluster number are valid.
Loops through all current clusters and loops through all points of the current cluster and sets the according
Dirichlet value.
NOTE: Sets Dirichlet value for Fermi level, i.e. only to be used in connection with
CurrentProblemC = 'solve-for-fermi-1' (LAPACK), 'solve-for-fermi-2'
(conjugate gradient)
Current in 1D
In 1D, there is an additional option to calculate the current by directly
integrating the current equation (MODULE
mod_integrate_current1D).
(Other option:
current_setup1D. Same as
current_setup2D
and
current_setup3D.)
The alternative method, which directly integrates the current equation is done
in subroutine
integrate_current1D.
SUBROUTINE integrate_current1D(chargeC,i1,i2,fermiV,E_FL,E_FR, &
curL,curR,curV,bndC,optionC)
Solves transport equation by using the full density.
j = mobility * density * grad(EF)
One has only two quasi-Fermi levels (one for electrons, one for holes).
For calculation of densities Fermi levels in MODULE fermi_level are valid.
Note: You must set all of them (for quantized states) to the value in fermi_nV.
Potential in MODULE potentials phiV is valid.
chargeC = 'el' or 'hl'
i1 , i2 ... >=1 , <= n_1 , i2>i1
fermiV (1..lengvecV)
Note: Fermi supposed to be continuous
at multiple points
but this can be changed to
discontinuous Fermi levels
... Output: Fermi level between points overwritten,
remains untouched in other regions.
E_FL,E_FR (Input/Output) ... Fermi level in i1, i2
curV (Input/Output) ... Current (particle current,
i1 --> i2 :: positive)
(1..n_1-1)
curV(n) = current between n and n+1 on
physical grid
curR,curL (Input/Output) ... Current that flows in or out.
bndC ... 'E_l_E_r' ... left/right Fermi level as boundary condition
'E_l_cur' ...
left Fermi level and current curL as boundary condition
'E_r_cur' ...
right fermilevel and current curR as boundary condition
Note: This subroutine only takes one iteration.
EFL, EFR, Fermi levels in [J]
cur in A/m2
NOTE:
Up to now, initial value for current is set to 1.0d0.
It would make sense to use curV as input to extract a starting value.
This method integrates the equation:
d/dx exp(E_F/kT) = J / (mu*k*T*kappa) , kappa
= n exp(-E_F/kT)
(which is obtained from j = mu * n * grad E_F)
The integration is carried out from i1 to i2 (i1,i2
on physical grid).
The current is given on material grid points (i.e. between physical
points) and it is stored in the variable
curV (1..n_1-1).
Recombination: div j = d/dx j = -R -->
j_R - j_L = - int R dx
The current in position i can be split up: curV(i) = cur_0 + j_rec
with j_recR - j_recL = - int R dx
NOTE: cur_0 is the current
that enters the device from the left contact.
Thus we get for points i and i+1:
exp(E_F/kT)|i+1 - exp(E_F/kT)| i = int_i-i+1 J
/ (mu*k*T*kappa) dx =
= int_i_i+1 cur_0 / (mu*k*T*kappa) dx + int_i_i+1 j_rec/ (mu*k*T*kappa)
dx
The added values of int_i-i+1 cur_0/ (mu*k*T*kappa) dx are
stored in exp_fermiV. (exp_fermiV(i1) = exp(E_FL/kT)),
the added values of int_i_i+1 j_rec/ (mu*k*T*kappa) dx
in exp_fermi_recV. (exp_fermi_recV(i1) = 0).
At first, cur_0 is set to 1.0d0. Then the integration is
carried out.
Finally cur_0 is scaled in a way that exp_fermiV
reaches the value exp(E_FR/kT) at the right which is the boundary
condition.
Note: exp_fermiV(i1) = exp(E_FL/kT)
--> If one begins at point i and wants to rescale cur_0
such that exp_fermiV(end_point) = exp(E_FR/kT):
exp(E_FR/kT) -
exp(E_Fi/kT) = scale * int_i_end cur_0 / (mu*k*T*kappa) dx
+ int_i_end j_rec / (mu*k*T*kappa) dx
--> scale = [exp(E_FR/kT) - exp(E_Fi/kT) - int_i_end j_rec
/ (mu*k*T*kappa) dx ] / [int_i_end cur_0 / (mu*k*T*kappa) dx ] =
[exp(E_FR/kT) - exp(E_Fi/kT) - (exp_fermi_recV(end_point)-exp_fermi_recV(i))
/ exp_fermiV(end_point) - exp_fermiV(i))
Recombination/generation and mobilities are obtained via the subroutines in
MODULE
get_transport_on_guent.
A positive recombination rate means generation of electrons/holes.
This subroutine is called in the internal subroutine
calculate_new_fermi_int_current in subroutine
solve_current_problem1D.
Matrix-vector product
The subroutines that calculate the matrix-vector products are:
subroutine
current_matvec,
subroutine
psolve_jac_current contained in
MODULE
current_routines.
|