Van de Walle (PRB 39, 1871 (1989)) describes the
dependence on strain of the valence and conduction band edges (zincblende). In
this paper strain is decribed along [001] and [111] directions. In a general
case, when we have [a b c] growth direction you will find the equations below (will
have to be written down!).
I made some derivations in order to understand what the single band model for
holes is, when we have an arbitrary growth direction. The results may be
summarized as follows (I follow a paper of Chao and Chuang PRB 46, 4110):
- 1) In the arbitrary growth direction strain tensor eij can have 6
independent components.
- 2) Starting point is the Hamiltonian of Bir and Pikus, that takes into
account eij.
- 3) We shall neglect split-off band.
- 4) In the case of [001] growth direction, when strain tensor has diagonal
components only with exx=eyy<>ezz, we find that there are two bands, each
of them is twofold degenerated. Effective mass tensor is diagonal in the
crystal system and does not depend on strain.
- 5) But in the case of arbitrary growth direction the situation is more
complicated:
- There are two bands, as it was before, each of them is twofold degenerated.
- Effective mass tensor is not diagonal in the crystal system. Its elements
depend on strain.
Taking into account all above, I'd like to ask you: Can we use this model in
nextnano3?
It seems to me, that we need three modifications only:
- 1) Correct dependence of the band extrema on strain are implemented
already.
- 2) Dependence of the effective mass on strain - not very difficult.
- 3) Nondiagonal effective mass tensor - this is more difficult, because it
requires coordinate transformation of such terms as kxky, kykz and kxkz. But
such problems must be solved already in k.p. As for the discretization of
nondiagonal terms, I think that this is already done, because even if in
the crystal system we don't have them, they appear in the simulation system.
In the MODULE mod_vb_abs_energy band edges are calculated starting from
k.p Hamiltonian. So, band edges are correct for any orientation. In order to get
the band edge shift due to strain, we
diagonalize numerically at each grid point the strain dependent
k.p Hamiltonian for k=0.
k.p parameters in the database are used to determine band shifts due
to strain (MODULE mod_vb_abs_energy). Program computes the k.p and strain
dependent Hamiltonian for each grid point, diagonalizes it, and takes the
eigenvalues for the band edges. So, band edges for single band calculations are
k.p parameter dependent.
Last days I studied how 1-band Schrödinger is implemented in nextnano3. I
tried to understand how band edges are determined. As far as I understood, in
order to determine shifts due to strain, you diagonalise k.p Hamiltonian for
every grid point. Eigenvalues are taken for the band edges.
Despite of the fact that FUNCTION vb_abs_energy(grid_point,band_number)
has
detailed comments that have helped me a lot, there are some points that I'd like
to clarify.
- 1) What are the basis vectors that are defined by
CALL
get_basis_vectors(strain_tensor_cxyz,basis_vec,Spin_orbit_HamM)?
- 2) Can you please explain how do you choose energy of "heavy" hole and
"light" hole among the eigenvalues. I'd be very happy to know just a physical
algorithm of yours.
I diagonalize the k.p Hamiltonian including deformation potentials at
k=0 at every spatial grid point. The eigenvalues are the local valence
band edges. The problem now is to find out which eigenvalue belongs to which
band:
What I do is, I project the eigenstates onto the spin-orbit eigenstates named as
basis_vec. One problem though is the quantization axis of these states which is
the z axis by default. Since the strain is breaking the symmetry the states are
getting mixed. I tried to overcome this problem for some extent by trying to
find some symmetry axis of the strain tensor. Then, anticipating that this is
the best possible quantization axis, I rotate the
basis_vec such that the
quantization axis coincides with the symmetry axis of the strain. Then I project
the eigenstates from the k.p diagonalization onto the rotated
basis_vec and just make a majority vote: If a state is 60% heavy and 40 % light hole I
call it "heavy hole". I checked this procedure on some structures but mainly in
1D, so maybe you can tell me if it works properly in your case.
I'm studying your subroutine vb_abs_energy.
Let me ask a "technical" question.
In the subroutine get_basis_vectors you rotate eigen vectors of the
spin-orbit
coupling Hamiltonian. Right?
In principle, eigenvector is a spinor. You independently rotate the first 3
components of it and the last 3 components of it. Am I right?
It looks like (maybe I'm wrong) that you don't rotate spin part of the wave
function.
As for me this is a mistake. As far as I remember from Landau & Lifshitz, when
you rotate axis you also mix different spins.
I made a routine that can calculate effective mass and band shifts of holes
for zinc-blende. My idea is very similar to what Matthias has done, but there
are some changes.
There are two possible situations:
- a) There is no strain, or strain ellipsoid is a sphere.
- b) otherwise
The problem: Diagonalize strain dependent k.p Hamiltonian and distinguish heavy,
light and split-off hole states.
Solution: Distinguish them by projecting them on the states |j,m>
Quantization axis is chosen in such a way:
case a) z-axis coincides with k vector;
case b) quantization axis coinsides with the principal axis of
the strain tensor that corresponds to the maximum eigenvalue.
In situation a) we cannot talk about effective-mass tensor because of
warping. However, it is possible to define effective mass along any direction.
I hope, of course, that k.p calculation will work fine.
But now I'd like to plug this routine into the code and see what will be the
result.
For that I need routine that can discretise 1/mxy terms.
Let me inform you about a test that I performed in order to understand if
band shifts due to strain are calculated correctly.
The structure that I consider is a 1D quantum well GaAs-4K/InAs-4K/GaAs-4K.
In order to check the program, I wrote my own subroutine and compared the band
shifts.
I considered 2 growth direction: [0 0 1] and [3 1 1].
Let us compare band edges in InAs (inside the well, where we have nonzero
strain):
In the [001] direction our results are the same:
Matthias, hh,lh,s/o 0.286451116426106 8.643207818113535E-002 -0.568763735470029
Michael, hh,lh,s/o 0.286451116426106 8.643207818113536E-002 -0.568763735470029
But for [311] growth direction they are different:
Matthias, hh,lh,s/o 0.304729495783759 9.995444713679591E-002
9.995444713679591E-002
Michael, hh,lh,s/o 9.995444713679575E-002
0.304729495783759 -0.561455237788657
Namely, your hh is equal to my lh, your lh and s/o are the same (!) and equal to
my hh. I'd like to say that your results are not accurate enough, on my opinion.
If you don't agree, please, let me know.
Basically, version of Matthias produced correct band edge for heavy and light
hole, but noncorrect for s/o band (it was the same as for the light hole).
Module that I sent you executes subroutine of Matthias, then executes my
subroutine, prints two values to compare and returns values of Matthias for the
function value.
So, it only compares.
I left it so because I didn't want to change code too much. You know that my
works for zincblende only yet.
When you have time, you change it in such a way that for wurzite it takes
Matthias value, for zb it takes mine.
As for me, I have my version where my values are taken instead of the values of
Matthias.
After your make a new release you should update tutorial 1D #5.
--
I send you my routine.
Actually, I wanted to spend some time and made also the modification for
wurtzite.
So, you will get:
- 1)
mod_kp2sb.f90
It contains subroutines that transfer parameters from your code to my routine,
subroutine that calculates band shifts and effective masses.
- 2) I send you
module vb_abs_energy.f90.
I added lines in order to compare your and my results:
! my routine
call init_strain(strain_tensor_cxyz) ! passes strain to my code
as = (2*ms1+ls1)/3 !
conversions to "traditional" parameters
bs = (ls1-ms1)/3
ds = ns1/sqrt(3d0)
print *,'d',d
call init_LMN_param( L,M,N,as,bs,ds,d) ! passes k.p - parameters to my code
call
get_shifts(hh_energy1,lh_energy1,so_energy1)
! calculates band edges
print
*,'-----------------------------------------------------'
print *,'Matthias, hh,lh,s/o',hh_energy, lh_energy, so_energy
print *,'Michael, hh,lh,s/o',hh_energy1,lh_energy1,so_energy1
print *,'-----------------------------------------------------'
! end of my routine
The second thing that I'd like to do is to implement effective mass
calculation.
Strain also enters into the single-band calculation.
There are two effects: Band shifts and effective mass changes.
As for the band shift I used my routine that is, presumably, correct. As for the
effective mass, I calculated their dependence on strain and put the values into
the input file. For that I needed to put them in the input file.
I believe that there should be a flag, that will allow the user to choose between
constant effective masses (taken from database) and effective masses calculated
from k.p parameters. Of course, one needs nondiagonal masses.
> So for 311 the strain is "more important" than for 001. Maybe this could be a reason for the difference.
> > I solved Schrödinger equation for a rather large well with 1D k.p. After that I calculated effective mass from k.p parameters and solved the same in the effective mass approximation. Since it was 1D, I could use isotropic masses.
> > The interesting thing was that for [001] direction results for the 1st state were almost the same, but for [311] they were very different.
> > May be it should be so.
> > I have an opinion, that when I touch something with [311] direction
> > it can be wrong. This was with strain and with band edges. We have to check
k.p.
main calls SUBROUTINE
gen_new_band_struc.
Lots of routines contain these expressions explicitly. By calling one of the
"never-called" functions, one should be able to "save" these lines in the code.