# [Pw_forum] Hubbard forces without Hubbard U

Sclauzero Gabriele gabriele.sclauzero at mat.ethz.ch
Fri Aug 15 11:37:10 CEST 2014

Dear Dmitry,

I am sure there is nothing puzzling in there. When lda_plus_u is true, but all Hubbard_U’s are zero, the Hubbard component of the forces are still computed and added to the total force, but such component will be zero (try verbosity=‘high’ to get them printed on output). In the subroutine force_hub you can find the line:

forceh(ipol,alpha) = forceh(ipol,alpha) -  v%ns(m2,m1,is,na) * dns(m1,m2,is,na)

which corresponds exactly to the formula you wrote here below. The point is that the Hubbard potential
v%ns  =\frac{U}{2} \sum_{Imm’} (\delta_{mm’} -2n^{I}_{m’m})
depends linearly on Hubbard_U, so all its components will be zero if U=0.

HTH

GS

> Dear Prof. Cococcioni,
>
> I use the QE 5.1 and It seems that the Hubbard forces are calculated well. But I'm a little puzzled by the next situation.
> An analytical expression for the Hubbard forces can be written as:
>
> F^{U}_{\alpha i} = - \frac{U}{2} \sum_{Imm'} (\delta_{mm'} -2n^{I}_{m'm}) \frac{\partial n^{I}_{mm'}}{\partial \tau_{\alpha i}}.
>
> Here alpha is a number of displacement atom in i-direction and U is a Hubbard correction.
> But I have not found the procedure of multiplication by the Hubbard U correction in the subroutines where Hubbard forces are calculated. Below the appropriate lines of the corresponding subroutines are listed:
>
> SUBROUTINE forces()
>   ...
>   CALL force_hub(forceh)
>   ...
>   force(ipol,na) = force(ipol,na) + forcenl(ipol,na) + forceion(ipol,na) + forcelc(ipol,na) + forcecc(ipol,na) + forcescc(ipol,na) + forceh(ipol,na)
>   ...
> END SUBROUTINE forces
>
>
> SUBROUTINE force_hub(forceh)
>   ...
>   forceh(ipol,alpha) = forceh(ipol,alpha) - v%ns(m2,m1,is,na) * dns(m1,m2,is,na)
>   ...
> END SUBROUTINE force_hub
>
> --
>
> Best regards,
> Dr. Dmitry Novoselov
>
> Institute for Metal Physics,
> Yekaterinburg, Russia
>
>
> 2014-08-14 19:43 GMT+06:00 Matteo Cococcioni <matteo at umn.edu>:
>
> Dear Dmitry,
>
> it should be working. what version of the code are you using?
>
> Matteo
>
>
> On Wed, Aug 13, 2014 at 10:47 AM, Dmitry Novoselov <dnovoselov at gmail.com> wrote:
> Dear all,
>
> It seems that in the part of the code where the Hubbard forces is calculated (forces and force_hub subroutines), the Hubbard U correction does not taken into account. Is it ok?
>
> Thank you!
>
> --
>
> Best regards,
> Dr. Dmitry Novoselov
>
> Institute for Metal Physics,
> Yekaterinburg, Russia
>
>
> _______________________________________________
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://pwscf.org/mailman/listinfo/pw_forum
>
>
> _______________________________________________
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://pwscf.org/mailman/listinfo/pw_forum
>
>
>
> --
>
> Best regards,
> Dr. Dmitry Novoselov
>
> Institute for Metal Physics,
> Yekaterinburg, Russia
>
> _______________________________________________
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://pwscf.org/mailman/listinfo/pw_forum

Dr. Gabriele Sclauzero
Materials Theory - ETHZ
ETH Zurich, HIT G 43.2
Wolfgang-Pauli-Str. 27
8093 Zürich, Switzerland

Phone +41 44 633 94 10
Fax +41 44 633 14 59
gabriele.sclauzero at mat.ethz.ch
www.theory.mat.ethz.ch