|
| |
Current-Poisson problem
The current-Poisson problem comprises the self-consistent solution for the
potential and the quasi-Fermi levels. This problem is treated block-iteratively.
The iteration consists of two main steps:
- Calculate the new quasi-Fermi level for the
given potential.
- Calculate the new potential for the
given quasi-Fermi levels (Schrödinger-Poisson)
This is done in subroutine
current_poisson.
For a given potential the new quasi-Fermi level is calculated with subroutine
solve_current_problem.
In 1D, it is done by default by directly integrating the current equation (subroutine
calculate_new_fermi_int_current). See documentation of
setup current problem for further information. In
2D/3D this is done by solving the current equation as described in the script.
Coupling of current and Poisson equation
As described in the script the potential and the quasi-Fermi levels can be
coupled. This is done in subroutine
current_poisson_coupled.
The Newton method calculates the root of the functional
calulate_functional_coupled.
One additionally has to provide the subroutines
newton_corr_coupled and
calc_grad_coupled.
These are contained in MODULE
functionals_newton_coupled.
!----------------------------------------------
SUBROUTINE calculate_functional_coupled(xV,fV)
!----------------------------------------------
Calculates the functional fV(xV) for the Newton method.
Input:
xV (1..dim_problem)
Output:
fV (1..dim_problem)
Note: The dimension of the total equation is dim_problem.
dim_problem = 3 <--> poissonL,fermi_nL,fermi_pL = .TRUE.
= 2 <--> two of them are .TRUE.
= 1 <--> one is .TRUE.
= 0 <-->
all .FALSE.
poissonL, fermi_nL, fermi_pL are global in this MODULE.
----------------------
The functional is defined as follows:
div D = -
div grad phi = rho
- div j_n = - div mu |n| grad E_Fn =
-R
- div j_p = - div mu |p| grad E_Fp = R
/ - div D - rho \
--> F[phi,E_Fn,E_Fp] = | - div j_n + R |
\ - div j_p - R /
!-----------------------------------------
SUBROUTINE newton_corr_coupled(dxV,xV,fV)
!-----------------------------------------
Given xV,fV
(input) this subroutine calculates one Newton correction step
dxV (output):
fdjac * dxV = - fV
dxV = -
(fdjac)-1 * fV
fdjac = DF =
∂F/∂x
= Jacobian matrix
Comment:
These subroutines are used for subroutine
newton.
calc_grad_coupled is always called directly after newton_corr_coupled.
Hence the Jacobian is not deallocated in newton_corr_coupled but in calc_grad_coupled.
The subroutine sets up the Jacobian by calling subroutine
define_jacobian_coupled.
The Jacobian can be accessed by the subroutines
matvec_DF,
Tmatvec_DF (transpose
of Jacobian),
psolve_jac_current_DF (Jabobi preconditioner). The Jacobian consists of
the following parts (see also script): [subroutine
define_jacobian_coupled]
The Jacobian has following structure:
phi
fermi_n
fermi_p
phi 1 2 3
2
3
fermi_n 4 10
4 5 8
9
fermi_p 6 13
12
6 7 11
1 --> Poisson matrix
2 --> d rho/d E_Fn = d n/d phi
3 --> d rho/d E_Fp = d p/d phi
4 --> 'diffEF' for electrons
5 --> 'matrix' for electrons
6 --> current matrix with mode 'diffEF' for holes
7 --> 'matrix' for holes
8 --> d_rec / d_E_Fn [drec_dEFn8V] div(j_n)
equation
9 --> d_rec / d_E_Fp [drec_dEFp9V] div(j_n)
equation
10 --> d_rec / d_phi [drec_dphi10V] div(j_n) equation
12 --> d_rec / d_E_Fn [drec_dEFn12V] div(j_p) equation
11 --> d_rec / d_E_Fp [drec_dEFp11V] div(j_p)
equation
13 --> d_rec / d_phi [drec_dphi13V] div(j_p)
equation
It should be: 8 = -12 , 9 = -11 , 10 = -13
The setup for the parts of the Jacobian is done as described above:
- Poisson matrix: Has already been set up at the beginning. Does not change
during the whole calculation.
d rho/d E_Fn = d n/d phi: is calculated using the subroutines
deriv_density and deriv_donor_density.d rho/d E_Fp = d p/d phi: is calculated using the subroutines
deriv_density and deriv_acceptor_density.- Current matrix with mode '
diffEF' for electrons / with mode 'matrix' for
electrons: In the linear current matrix the density is considered as
coefficients. For calculation of the derivative one has to take into account
that the density itself is a function of the potential/Fermi level. So one has
to apply the chain rule. One thus gets two components for the derivative as
described in the script:
d/d_E_Fn (M.E_Fn)=M+D_E_FnM.
With mode 'matrix' M is defined, with mode 'diffEF' D_E_FnM.
See below (documentation of subroutine currentmatrix_for_coupled
for
further information). These four matrices are stored in MODULE current1D
(Note: They are not stored as symmetric matrices).
- Recombination is defined on each grid point and it only depends on
potential and Fermi level on this grid point. Hence its Jacobian is a simple
diagonal matrix.
!---------------------------------------------
SUBROUTINE calc_grad_coupled(xV,fV,gradientV)
!---------------------------------------------
calculates fVT * Df = gradientV
xV ... DIMENSION(1:lengvecV)
fV ... DIMENSION(1:lengvecV)
gradientV ... DIMENSION(1:lengvecV)
Comment:
These subroutines are used for subroutine
newton.
calc_grad_coupled is always called directly after newton_corr_coupled.
Hence the Jacobian is not deallocated in newton_corr_coupled but in calc_grad_coupled.
The Jacobian has been set up before: It can be accessed via subroutine
Tmatvec_DF (Transpose of Jacobian).
Setup of matrices for coupled solution
The routines for the setup of the four block matrices for the current equation
(see above) is done by the subroutines currentmatrix_for_coupled,
DIRICHLET_coupled1D, set_cur_differential1D and set_dirichlet_current(coupledL=.TRUE.) in
MODULE
mod_set_dirichlet_current.
See script for further information. The way of the setup is in principle the same
as for the linear current equation.
Predictor density
The functional is solved with a predictor density, so the wavefunctions are
considered fixed during the minimization of the functional.
In order to calculate the shift of the eigenvalues, the potential for which
the eigenvalues have been solved is stored in potV_old in MODULE
pass_over_mod.
|
