Tutorial 4: Polarization in GaN

 

Here we compute the polarization difference between wurtzite (W) and zinc-blende (ZB) structure of GaN. This polarization difference is responsible for the presence of a build-in electric field in heterostructures of group-III nitrides. The two structures (wurtzite and zinc-blende) are shown below.

Our basic steps will include obtaining SCF electron density for these structures followed by calculation of polarization and their subsequent difference

Ps = P(ZB) - P(W)

Note: Ps is a property of the interface rather than GaN alone.

Wurtzite GaN

Let's create a work directory GaN-W

1. Creating the structure file GaN-W.struct

GaN 4 atoms
H   LATTICE,NONEQUIV.ATOMS:  2186_P63mc
MODE OF CALC=RELA unit=bohr
  5.963131  5.963131  9.722374 90.000000 90.000000120.000000
ATOM  -1: X=0.33333333 Y=0.66666666 Z=0.50000000
          MULT= 2          ISPLIT= 4
      -1: X=0.66666667 Y=0.33333333 Z=0.00000000
Ga1        NPT=  781  R0=0.00005000 RMT=    1.9300   Z: 31.0
LOCAL ROT MATRIX:    1.0000000 0.0000000 0.0000000
                     0.0000000 1.0000000 0.0000000
                     0.0000000 0.0000000 1.0000000
ATOM  -2: X=0.33333333 Y=0.66666666 Z=0.87639300
          MULT= 2          ISPLIT= 4
      -2: X=0.66666667 Y=0.33333333 Z=0.37639300
N 1        NPT=  781  R0=0.00010000 RMT=    1.5000   Z:  7.0
LOCAL ROT MATRIX:    1.0000000 0.0000000 0.0000000
                     0.0000000 1.0000000 0.0000000
                     0.0000000 0.0000000 1.0000000
  12      NUMBER OF SYMMETRY OPERATIONS
 1 0 0 0.00000000
 0 1 0 0.00000000
 0 0 1 0.00000000
       1
...

It is the best to get the file from a repository BerryPI/tutorials/tutorial4/GaN-W/GaN-W.struct. Here the lattice parameters as well as atomic positions are fully optimized.

2. Initialize calculation

init_lapw -b -vxc 5 -rkmax 7 -numk 300

Here we opt for LDA exchange-correlation functional and the resultant 8x8x4 k-mesh, which is dense enough for the purpose of this calculations.

3. Run SCF calculation

run_lapw

It might be a good idea to speed up the process by getting more cores involved export OMP_NUM_THREADS=4

4. Run BerryPI

berrypi -k 8 8 8

The results of the calculation are

...
SUMMARY OF POLARIZATION CALCULATION
=======================================================================================
Value                           |  spin   |    dir(1)    |    dir(2)    |    dir(3)
---------------------------------------------------------------------------------------
Electronic polarization (C/m2)     sp(1)  [ 2.223396e-10,  1.205149e-07, -5.612307e-02]
Ionic polarization (C/m2)          sp(1)  [ 1.750356e-10, -2.049427e-07, -4.386009e-01]
Tot. spin polariz.=Pion+Pel (C/m2) sp(1)  [ 3.973752e-10, -8.442784e-08, -4.947239e-01]
---------------------------------------------------------------------------------------
TOTAL POLARIZATION (C/m2)          both   [ 3.973752e-10, -8.442784e-08, -4.947239e-01]
=======================================================================================
...

As we can see the polarization (-0.4947 C/m2) is non-zero only along Z-axis, which corresponds to [0001] direction.

Zinc-blende GaN

Now we can switch to another work directory GaN-ZB

1. Creating the structure file GaN-ZB.struct

GaN cubic
H   LATTICE,NONEQUIV.ATOMS:  6156_P3m1
MODE OF CALC=RELA unit=bohr
  5.963131  5.963131 14.606628 90.000000 90.000000120.000000
ATOM  -1: X=0.00000000 Y=0.00000000 Z=0.00000000
          MULT= 1          ISPLIT= 4
Ga1        NPT=  781  R0=0.00005000 RMT=    1.9300   Z: 31.0
LOCAL ROT MATRIX:    1.0000000 0.0000000 0.0000000
                     0.0000000 1.0000000 0.0000000
                     0.0000000 0.0000000 1.0000000
ATOM  -2: X=0.00000000 Y=0.00000000 Z=0.25000000
          MULT= 1          ISPLIT= 4
N 1        NPT=  781  R0=0.00010000 RMT=    1.5000   Z:  7.0
LOCAL ROT MATRIX:    1.0000000 0.0000000 0.0000000
                     0.0000000 1.0000000 0.0000000
                     0.0000000 0.0000000 1.0000000
ATOM  -3: X=0.33333333 Y=0.66666667 Z=0.33333333
          MULT= 1          ISPLIT= 4
Ga2        NPT=  781  R0=0.00005000 RMT=    1.9300   Z: 31.0
LOCAL ROT MATRIX:    1.0000000 0.0000000 0.0000000
                     0.0000000 1.0000000 0.0000000
                     0.0000000 0.0000000 1.0000000
ATOM  -4: X=0.33333333 Y=0.66666667 Z=0.58333333
          MULT= 1          ISPLIT= 4
N 2        NPT=  781  R0=0.00010000 RMT=    1.5000   Z:  7.0
LOCAL ROT MATRIX:    1.0000000 0.0000000 0.0000000
                     0.0000000 1.0000000 0.0000000
                     0.0000000 0.0000000 1.0000000
ATOM  -5: X=0.66666667 Y=0.33333333 Z=0.66666667
          MULT= 1          ISPLIT= 4
Ga3        NPT=  781  R0=0.00005000 RMT=    1.9300   Z: 31.0
LOCAL ROT MATRIX:    1.0000000 0.0000000 0.0000000
                     0.0000000 1.0000000 0.0000000
                     0.0000000 0.0000000 1.0000000
ATOM  -6: X=0.66666667 Y=0.33333333 Z=0.91666667
          MULT= 1          ISPLIT= 4
N 3        NPT=  781  R0=0.00010000 RMT=    1.5000   Z:  7.0
LOCAL ROT MATRIX:    1.0000000 0.0000000 0.0000000
                     0.0000000 1.0000000 0.0000000
                     0.0000000 0.0000000 1.0000000
   6      NUMBER OF SYMMETRY OPERATIONS
 1 0 0 0.00000000
 0 1 0 0.00000000
 0 0 1 0.00000000
       1
...

It is the best to get the file from a repository BerryPI/tutorials/tutorial4/GaN-ZB/GaN-ZB.struct. The structure has 6 atoms vs 4 atoms for the wurtzite structure. Here the lattice parameter a = b = 5.963131 Bohr is the same as for the wurtzite structure above (otherwise how do we imagine the heterointerface), however we cannot use the same c, even re-scaled. The ideal zinc-blende c = a*sqrt(6) = 14.606628 Bohr. Next steps are identical to the wirtzite case.

2. Initialize calculation

init_lapw -b -vxc 5 -rkmax 7 -numk 200

Here we opt for 200 k-points that results in 8x8x2 k-mesh. This is a consequence of smaller Brillouin zone for the ZB structure in comparison to W structure.

3. Run SCF calculation

run_lapw

4. Run BerryPI

berrypi -k 8 8 8

Here we keep k-mesh used for the Berry phase integration identical to the wurtzite case. The results of the calculation are

...
SUMMARY OF POLARIZATION CALCULATION
=======================================================================================
Value                           |  spin   |    dir(1)    |    dir(2)    |    dir(3)
---------------------------------------------------------------------------------------
Electronic polarization (C/m2)     sp(1)  [-1.219101e-03,  1.219101e-03, -1.727041e-04]
Ionic polarization (C/m2)          sp(1)  [ 1.165061e-10,  1.165088e-10, -4.644818e-01]
Tot. spin polariz.=Pion+Pel (C/m2) sp(1)  [-1.219101e-03,  1.219101e-03, -4.646545e-01]
---------------------------------------------------------------------------------------
TOTAL POLARIZATION (C/m2)          both   [-1.219101e-03,  1.219101e-03, -4.646545e-01]
=======================================================================================
...

This result can be quite puzzling, since we would expect the polarization to vanish in the ZB structure for symmetry reasons. The explanation is provided by David Vanderbilt and R. D. King-Smith Phys. Rev. B 48, 4442 (1993), Sec. III.E.

Polarization difference

Finally, we can compute the polarization difference

Ps = -4.646545e-01 - (-4.947239e-01) = 0.030 C/m2

The obtained result agrees with Ps = -0.029 C/m2 reported by Bernardini et al. Phys. Rev. B 56, R10024 (1997); it is also consistent with Ps = -0.022 C/m2 measured by Lähnemann et al. Phys. Rev. B 86, 081302(R) (2012) for the heterointerface W/ZB GaN. The sign of P and Ps depends on the orientation of our structure and will flip upon flipping the structure along c-axis.

Note 1: completed using WIEN2k 13.1, WIEN2WANNIER 1.0-beta2, BerryPI 1.2, Python 2.7.4, Numpy 1.6.2

Note 2: checked with WIEN2k 18.2, WIEN2WANNIER 2.0, BerryPI 1.4, Python 2.7.14, Numpy 1.16.0

Note 3: Does not work with Python 3.*

 

 https://github.com/spichardo/BerryPI/tree/master/tutorials/tutorial4