Possible symmetry error in full-DFPT nonlinear optical tensor for wurtzite ZnO

Possible symmetry error in full-DFPT nonlinear optical tensor for wurtzite ZnO, ABINIT 10.6.5

I am using ABINIT 10.6.5 to calculate Raman tensors and the electro-optic / nonlinear optical correction for wurtzite ZnO. I am using norm-conserving pseudo-potentials. I had an issue with ixc=11, as this is the exchange type used by the pseudo-potentials. I set ixc to 7, to allow the code to pass through the optdriver=5 checks. I have seen the same issue with PAW pseudo potentials.

The calculation uses the full-DFPT route for third derivatives, not PEAD:

optdriver7        5
usepead7          0
d3e_pert1_elfd7   1
d3e_pert1_phon7   1
d3e_pert1_atpol7  1 4
d3e_pert2_elfd7   1
d3e_pert3_elfd7   1

The workflow is:

Dataset 1: ground-state SCF
Dataset 2: NSCF wavefunctions, kptopt=2
Dataset 3: d/dk wavefunctions
Dataset 4: phonon + electric-field + strain first-order responses
Dataset 5: d/dkdk second-order wavefunctions
Dataset 6: d/dkde second-order wavefunctions
Dataset 7: third derivatives, optdriver=5

I then merged the second- and third-derivative DDBs:

ramano_DS4_DDB
ramano_DS7_DDB

to produce:

analysis.ddb

and analysed it with anaddb using:

dieflag    1
nlflag     1
elaflag    3
piezoflag  3
instrflag  1

ramansr    1
alphon     1
prtmbm     1

ABINIT correctly identifies the structure as wurtzite ZnO:

space group P6_3 m c (#186)
nsym = 12

However, the nonlinear optical tensor printed in both raman.abo and analysis.abo does not have the symmetry required by point group 6mm.

The tensor printed by anaddb is:

Non-linear optical coefficients d (pm/V)

   -0.000000   -0.000000   -0.000000  -22.617400   92.228497   -0.000000
   -0.000000   -0.000000   -0.000000  222.810118  -22.617400   -0.000000
   92.228497  222.810118  -37.956596   -0.000000   -0.000000  -22.617400

Assuming the usual contracted order:

xx yy zz yz xz xy

this gives, for example:

d_xxz =  92.228497 pm/V
d_yyz = 222.810118 pm/V
d_zxx =  92.228497 pm/V
d_zyy = 222.810118 pm/V
d_xyz = d_xzy = d_yxz = d_yzx = d_zxy = d_zyx = -22.617400 pm/V

For wurtzite ZnO, space group P6_3mc, point group 6mm, with z along the hexagonal c axis, the tensor should have the form:

        xx    yy    zz    yz    xz    xy
x        0     0     0     0   d15     0
y        0     0     0   d15    0      0
z      d31   d31   d33    0     0      0

There should be 3 unique non-zero values, but in the results above there are 4. Also:

d_xxz should equal d_yyz
d_zxx should equal d_zyy
d_xyz, d_xzy, d_yxz, d_yzx, d_zxy, d_zyx should vanish

The reported tensor violates these conditions by tens to more than 100 pm/V, so this does not look like ordinary numerical noise.

There is one convergence warning in the ABINIT output:

nstep=200 was not enough SCF cycles to converge;
maximum residual=2.141E-22 exceeds tolwfr=1.000E-22

but this is only a factor of about two above an extremely tight tolerance and does not seem sufficient to explain such a large symmetry violation. It is however, strange that convergence at this steps is apparently instantaneous in one step and after the first step there is no improvement in the result. If this is expected then the number of steps could be reduced for this step.

Could this indicate a problem in the symmetrisation, Cartesian transformation, storage in the DDB, or anaddb reconstruction of the third-order electric-field tensor for full-DFPT optdriver=5, usepead=0 calculations? Or have I made a mistake with the input. Or maybe I am interpreting the output incorrectly. I would be very grateful for any help.
Thanks
John

I can provide the following files:

raman.abi
raman.abo
ramano_DS4_DDB
ramano_DS7_DDB
analysis.ddb
analysis.abi
analysis.abo
analysis.files

analysis.abo (43.9 KB)

raman.abi (6.7 KB)

raman.abo (319.7 KB)

analysis.abi (1.3 KB)

Hello John,

My very first question is : do you see the same problem of symmetry with the PEAD method? If you have tried it would be useful to know.

Concerning ixc, nonlinear computations with GGA functionals are not implemented, so only LDA is authorized by Abinit.

About the convergence warning : it seems that in your case the numerical noise does not allow a convergence of the 2nd WFs derivatives with tolwfr=10E-22. It happens, depending on the system and the convergence with respect to cutoff and kpoints (highly converged computations require bigger WFs which are more difficult to converge). At this step the Sternheimer equation is non self-consistent so if the maximal precision of the wavefunctions is reached at the first step, the other iterations won’t change anything, as you can see. To avoid these useless iterations increase tolwfr for dataset 5 and 6 (10E-21 or 10E-20 should be fine). I agree with you that this level of precision should be more than sufficient to give correct 3rd order derivatives.

Finally, there is a on-going investigation on the correctness of the DDB produced by the nonlinear driver (optdriver=5). The cartesian tensors written in the abinit output file are correct, but maybe the tensor in reduced coordinates in the DDB file may be wrong. An other open question (for now) is if anaddb compensates the possible error in the DDB files.

If there is an error, it will be corrected in the next release of the production version (so probably in the future 10.8.3).

Best,
Lucas Baguet

Thanks for the reply.

Yes, I did try the PEAD method and it gave me the symmetry I was expecting.

For the time being I will stick with the PEAD approach and look again when a new release comes out.

All the best
John