SCF calculation with SOC fails to converge.

[jingan-181]

Describe the bug

Using ABACUS v3.2.0 to perform the self-consistent calculation of monolayer SnTe with SOC, the charge density failed converge. The cause of the problem is speculated that the total magnetic moment is not 0 when the electron iterates. I set all atomic magnetic moments to 0 in the STRU file. The SnTe should be a non-magnetic material.
SnTe-scf.zip

Expected behavior

No response

To Reproduce

No response

Environment

No response

Additional Context

No response

[Satinelamp]

Just out of curiosity, I did some tests.
Actually, if you do more electronic step iterations, the convergence of DRHO would get stuck to about 10-7, which might be caused by the numerical error of EDIFF, as the log shows after the error of DRHO reach to 10-7, the total energy oscillates at the scale of 10-12. But usually, when the error of total energy reaches 10-7, the error of DRHO also reaches 10-7. This failure might be related to the numerical error caused by the precision of floating point number, which decides that the error of total energy couldn't be lower.

START CHARGE      : atomic
 DONE(11.7881    SEC) : INIT SCF
 ITER   ETOT(eV)       EDIFF(eV)      DRHO       TIME(s)    
 GE1    -9.991340e+03  0.000000e+00   2.934e-02  5.206e+00  
 GE2    -9.991425e+03  -8.482276e-02  1.101e-02  5.138e+00  
 GE3    -9.991420e+03  4.438911e-03   4.542e-03  5.085e+00  
 GE4    -9.991417e+03  3.527439e-03   2.514e-03  5.079e+00  
 GE5    -9.991414e+03  2.498785e-03   1.744e-03  5.101e+00  
 GE6    -9.991414e+03  3.131484e-04   1.136e-03  5.132e+00  
 GE7    -9.991414e+03  -1.932493e-05  8.682e-04  5.100e+00  
 GE8    -9.991414e+03  -6.721802e-05  7.221e-04  5.090e+00  
 GE9    -9.991414e+03  -1.373189e-04  6.538e-04  5.123e+00  
 GE10   -9.991414e+03  -1.333528e-04  5.884e-04  5.084e+00  
..........................
 GE1980 -9.991415e+03  1.113688e-10   4.777e-07  5.081e+00  
 GE1981 -9.991415e+03  1.546789e-12   3.736e-07  5.077e+00  
 GE1982 -9.991415e+03  5.444697e-10   4.739e-07  5.074e+00  
 GE1983 -9.991415e+03  -3.093578e-12  5.017e-07  5.077e+00  
 GE1984 -9.991415e+03  -3.712293e-11  5.461e-07  5.081e+00  
 GE1985 -9.991415e+03  -4.547559e-10  4.878e-07  5.077e+00  
 GE1986 -9.991415e+03  5.645779e-10   5.147e-07  5.080e+00  
 GE1987 -9.991415e+03  -4.563027e-10  4.475e-07  5.085e+00  
 GE1988 -9.991415e+03  1.036348e-10   4.992e-07  5.086e+00  
 GE1989 -9.991415e+03  4.021651e-11   5.128e-07  5.077e+00  
 GE1990 -9.991415e+03  -2.877027e-10  4.021e-07  5.084e+00  
 GE1991 -9.991415e+03  9.435412e-11   3.404e-07  5.114e+00  
 GE1992 -9.991415e+03  2.026293e-10   5.650e-07  5.094e+00  
 GE1993 -9.991415e+03  -3.588550e-10  3.698e-07  5.083e+00  
 GE1994 -9.991415e+03  1.577725e-10   5.102e-07  5.079e+00  
 GE1995 -9.991415e+03  8.971375e-11   4.855e-07  5.079e+00  
 GE1996 -9.991415e+03  2.768752e-10   4.869e-07  5.082e+00  
 GE1997 -9.991415e+03  -3.573082e-10  4.130e-07  5.118e+00  
 GE1998 -9.991415e+03  3.201853e-10   5.453e-07  5.101e+00  
 GE1999 -9.991415e+03  -8.197981e-11  5.451e-07  5.120e+00  
 GE2000 -9.991415e+03  -1.670532e-10  4.134e-07  5.106e+00  
 !! CONVERGENCE HAS NOT BEEN ACHIEVED !!

Here is the log if I run 2000 iterations.
log.txt

Here is another successful relax type of calculation.

 STEP OF RELAXATION : 74
 -------------------------------------------
 DONE(1.301e+03  SEC) : INIT SCF
 ITER   ETOT(eV)       EDIFF(eV)      DRHO       TIME(s)    
 GE1    -2.239390e+04  0.000000e+00   3.875e-04  1.174e+00  
 GE2    -2.239388e+04  1.832054e-02   9.978e-03  1.235e+00  
 GE3    -2.239390e+04  -1.857396e-02  2.176e-03  1.235e+00  
 GE4    -2.239391e+04  -1.892934e-03  1.152e-04  1.236e+00  
 GE5    -2.239391e+04  -8.096508e-06  4.857e-05  1.236e+00  
 GE6    -2.239391e+04  -5.891316e-07  4.472e-05  1.236e+00  
 GE7    -2.239391e+04  2.479162e-07   3.807e-05  1.235e+00  
 GE8    -2.239391e+04  2.609154e-07   2.487e-05  1.239e+00  
 GE9    -2.239391e+04  6.691130e-07   7.984e-06  1.115e+00  
 LARGEST GRAD (eV/A)  : 1.641e-02
 -------------------------------------------
 STEP OF RELAXATION : 75
 -------------------------------------------
 DONE(1.315e+03  SEC) : INIT SCF
 ITER   ETOT(eV)       EDIFF(eV)      DRHO       TIME(s)    
 GE1    -2.239390e+04  0.000000e+00   2.420e-04  1.174e+00  
 GE2    -2.239390e+04  4.234243e-03   3.071e-03  1.237e+00  
 GE3    -2.239390e+04  -3.436078e-03  1.868e-03  1.237e+00  
 GE4    -2.239391e+04  -1.362326e-03  7.110e-05  1.235e+00  
 GE5    -2.239391e+04  -3.425453e-06  2.836e-05  1.235e+00  
 GE6    -2.239391e+04  -2.739889e-07  2.937e-05  1.237e+00  
 GE7    -2.239391e+04  2.339735e-07   2.278e-05  1.235e+00  
 GE8    -2.239391e+04  1.194795e-06   9.541e-06  1.113e+00  
 LARGEST GRAD (eV/A)  : 8.838e-03

[Satinelamp]

the test in #1945 shows the difference of etot between different repeated simulations is about 10^-10.

[jingan-181]

Using HSE+SOC calculation will not encounter this problem.

HSE-scf.zip

[dyzheng]

I have checked this case of SnTe, HR of nonlocal part is not Hermitian for local matrix, which will lead to non-Hermitian for HK. This is not expected behavier, but coefficient of effective D is right with formula. I will check the case of HSE to double check this problem.

[wenfei-li]

I'll have a look into this issue

[wenfei-li]

first of all, I checked and the hamiltonian matrix is hermitian; I also wrote a separate script for calculating the soc vnl matrix, and the result matches that from ABACUS; finally, I tried checking out previous versions of ABACUS, and the magnetic moment has always been present
so it seems to me that this is not a bug, rather it is because ABACUS has poor convergence algorithm for nspin = 2 and 4

[wenfei-li]

btw the code for preparing the nonlocal part of soc pseudopotential in ABACUS is rather difficult to understand imo (lack of comment, strange way of arranging the indexes, etc.), I attach my script here. It does the equivalent of generating f_coeff in ABACUS.
soc.zip