scnc
/
quantum-espresso
mirror of https://gitlab.com/QEF/q-e.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Ford-PHonon - part 29 - fixes

This commit is contained in:
fabrizio22 2021-09-02 11:52:17 +02:00
parent 04b5480a8d
commit 4d7a25a182
12 changed files with 76 additions and 64 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

14
LR_Modules/adddvhubscf.f90
View File

@ -12,14 +12,18 @@ SUBROUTINE adddvhubscf (ipert, ik)
!
!! DFPT+U
!! This routine calculates the SCF derivative of the Hubbard potential times \(\psi\):
!! $$ |\Delta V_{SCF}(k+q,is) \psi(\text{ibnd},k,is)\rangle =
!! - \sum_{I,m1,m2} \text{Hubbard}_U(I)\cdot \text{dnsscf}(m1,m2,is,I,\text{imode})\cdot
!! |S\phi(I,k+q,m1)\rangle\langle S\phi(I,k,m2|\psi(\text{ibnd},k,is)\rangle $$
!! \begin{equation}\notag
!! \begin{split}
!! |\Delta V_{SCF}(k+q,is) \psi(\text{ibnd},k,is)\rangle =
!! - &\sum_{I,m1,m2} \text{Hubbard_U}(I)\cdot\text{dnsscf}(m1,m2,is,I,\text{imode})\cdot \\
!! &|S\phi(I,k+q,m1)\rangle\langle S\phi(I,k,m2|\psi(\text{ibnd},k,is)\rangle
!! \end{split}
!! \end{equation}
!
!! Addition of the \(text{J0}\) terms:
!! Addition of the \(\text{J0}\) terms:
!
!! $$ + \sum_{I,m1,m2} \text{Hubbard_J0}(I)\cdot \text{dnsscf}(m1,m2,\text{isi},I,\text{imode})\cdot
!! \|S\phi(I,k+q,m1)\rangle\langle S\phi(I,k,m2)\|\psi(\text{ibnd},k,is)\rangle $$
!! |S\phi(I,k+q,m1)\rangle\langle S\phi(I,k,m2)|\psi(\text{ibnd},k,is)\rangle $$
!
!! Where:
!! \(\text{is}\) = current_spin;

4
PHonon/PH/alpha2f.f90
View File

@ -140,7 +140,7 @@ END SUBROUTINE read_polarization
!--------------------------------------------------------------------
SUBROUTINE read_lam()
!------------------------------------------------------------------
!! This routine reads \(\text{lambda}_{q nu}\) & \(\text{omega}_{q nu}\)
!! This routine reads \(\text{lambda}_\text{q nu}\) & \(\text{omega}_\text{q nu}\)
!! from lambda*.dat
!
USE kinds, ONLY : DP
@ -296,7 +296,7 @@ END SUBROUTINE compute_a2F
!-----------------------------------------------------------------
SUBROUTINE compute_lambda()
!---------------------------------------------------------------
!! This routine computes \(\text{omega_ln}\) & \(\text{lambda}\).
!! This routine computes \(\text{omega_{ln}}\) & \(\text{lambda}\).
!
USE kinds, ONLY : DP
USE modes, ONLY : nmodes

10
PHonon/PH/delta_sphi.f90
View File

@ -12,7 +12,7 @@ SUBROUTINE delta_sphi (ikk, ikq, na, icart, nah, ihubst, wfcatomk_, wfcatomkpq_,
dqsphi, dmqsphi, iflag)
!---------------------------------------------------------------------------------
!
! DFPT+U: This routine calculates a vector at k:
!! DFPT+U: This routine calculates a vector at k
!
! |\Delta_{-q}(S_{k+q} \phi_(k+q,I,m)) > = S_{k} | \d^{icart}phi_(k,nah,m) > +
! \sum{l1,l2} [ | \dbeta^{icar_}_(k,na_,l1) > * qq_nt(na_, l1 ,l2) *
@ -20,7 +20,7 @@ SUBROUTINE delta_sphi (ikk, ikq, na, icart, nah, ihubst, wfcatomk_, wfcatomkpq_,
! | \beta_(k,na_,l1)> * qq_nt(na_, l1 ,l2) *
! < \dbeta^{icar_}_(k+q,na_,l2) | phi_(k+q,nah,m) > ]
!
!! and also a vector at k+q :
!! and also a vector at k+q.
!
! |\Delta_q(S_{k} \phi_(k,I,m)) > = S_{k+q}| \d^{icart}phi_(k+q,nah,m) > +
! \sum{l1,l2} [ | \dbeta^{icar_}_(k+q,na_,l1) > * qq_nt(na_, l1 ,l2) *
@ -32,8 +32,10 @@ SUBROUTINE delta_sphi (ikk, ikq, na, icart, nah, ihubst, wfcatomk_, wfcatomkpq_,
! |\Delta_{-q}(S_{k+q} \phi_(k+q,I,m)) >
! iflag = 0 : calculate ONLY |\Delta_{-q}(S_{k+q} \phi_(k+q,I,m)) >
!
! Written by A. Floris
! Modified by I. Timrov (01.10.2018)
!! See source comment for details on the implemented formulas.
!
!! Written by A. Floris.
!! Modified by I. Timrov (01.10.2018).
!
USE kinds, ONLY : DP
USE uspp_param, ONLY : nh, nhm

5
PHonon/PH/dvqpsi_us_only.f90
View File

@ -9,8 +9,9 @@
!----------------------------------------------------------------------
subroutine dvqpsi_us_only (ik, uact, becp1, alphap)
!----------------------------------------------------------------------
!! This routine calculates dV_bare/dtau * psi for one perturbation
!! with a given q. The displacements are described by a vector uact.
!! This routine calculates \(\text{dV_bare}/\text{dtau}\cdot\text{psi}\)
!! for one perturbation with a given q.
!! The displacements are described by a vector uact.
!! The result is stored in \(\text{dvpsi}\). The routine is called for
!! each k-point and for each pattern u. It computes simultaneously all
!! the bands.

50
PHonon/PH/dvscf_interpolate.f90
View File

@ -10,8 +10,8 @@ MODULE dvscf_interpolate
!----------------------------------------------------------------------------
!!
!! Module for Fourier interpolation of phonon potential dvscf.
!! Uses the output of dvscf_q2r.x program to compute the induced part of
!! dvscf at each q point.
!! Uses the output of \(\texttt{dvscf_q2r.x}\) program to compute the induced
!! part of \(\text{dvscf}\) at each q point.
!
! See header of dvscf_q2r.f90 for details.
!
@ -61,8 +61,8 @@ MODULE dvscf_interpolate
!----------------------------------------------------------------------------
SUBROUTINE dvscf_interpol_setup()
!--------------------------------------------------------------------------
!! Do setups for dvscf_r2q. Called in phq_setup at the beginning of each
!! q point calculation.
!! Do setups for \(\texttt{dvscf_r2q}\). Called in \(\texttt{phq_setup}\)
!! at the beginning of each q-point calculation.
!--------------------------------------------------------------------------
!
USE kinds, ONLY : DP
@ -234,20 +234,20 @@ MODULE dvscf_interpolate
!----------------------------------------------------------------------------
SUBROUTINE dvscf_r2q(xq, u_in, dvscf)
!--------------------------------------------------------------------------
!!
!! Read inverse Fourier transformed potential w_pot (written by
!! dvscf_q2r.x) and Fourier transform to compute dvscf at given q point.
!!
!! Originally proposed by [1], long-range part described in [2].
!! [1] Eiguren and Ambrosch-Draxl, PRB 78, 045124 (2008)
!! \(\texttt{dvscf_q2r.x}\)) and Fourier transform to compute
!! \(\text{dvscf}\) at given q point.
!
!! Originally proposed by [1], long-range part described in [2].
!! [1] Eiguren and Ambrosch-Draxl, PRB 78, 045124 (2008)
!! [2] Xavier Gonze et al, Comput. Phys. Commun., 107042 (2019)
!!
!! dvscf(r,q) = exp(-iqr) (dvlong(r,q) + sum_R exp(iqR) w_pot(r,R))
!!
!
!! $$ \text{dvscf}(r,q) = \text{exp}(-iqr) (\text{dvlong}(r,q) +
!! \sum_R \text{exp}(iqR) \text{w_pot}(r,R)) $$
!
!! In this subroutine, pool parallelization is used to distribute R points
!! to nodes so that the root of each pool reads different w_pot file
!! simultaneously, reducing the io time.
!!
!--------------------------------------------------------------------------
!
USE kinds, ONLY : DP
@ -473,15 +473,14 @@ MODULE dvscf_interpolate
!----------------------------------------------------------------------------
SUBROUTINE dvscf_shift_center(dvscf_q, xq, shift_half, sign)
!----------------------------------------------------------------------------
!! Shift center of phonon potential.
!! If sign = +1, shift center from origin to tau (used in dvscf_r2q).
!! Shift center of phonon potential.
!! If sign = +1, shift center from origin to tau (used in dvscf_r2q).
!! If sign = -1, shift center from tau to origin (used in dvscf_q2r).
!!
!
!! For ipol = 1, 2, 3, if shift_half(ipol) is true, the origin is set at
!! r = 0.5 for direction ipol. If false, the origin is set at r = 0.0.
!! Setting shift_half = .true. is useful when the size of supercell is odd,
!! becuase the center of the supercell is at (0.5, 0.5, 0.5).
!!
!----------------------------------------------------------------------------
USE kinds, ONLY : DP
USE constants, ONLY : tpi
@ -573,19 +572,20 @@ MODULE dvscf_interpolate
!----------------------------------------------------------------------------
!! This subroutine calculates the long-range dipole potential for given
!! xq and zeu.
!! Input xq: q vector in Cartesian coordinate
!!
!! Input xq: q vector in Cartesian coordinate.
!
!! Currently, only the dipole part (Frohlich) is implemented. The quadrupole
!! potential is not implemented.
!!
!
!! [1] Xavier Gonze et al, Comput. Phys. Commun., 107042 (2019)
!!
!
!! Taken from Eq.(13) of Ref. [1]
!! dvlong(G,q)_{a,x} = 1j * 4pi / Omega * e^2
!! * [ (q+G)_y * Zstar_{a,yx} * exp(-i*(q+G)*tau_a)) ]
!! / [ (q+G)_y * epsilon_yz * (q+G)_z ]
!! $$ \text{dvlong}(G,q)_{a,x} = 1j \cdot 4\pi / \text{omega} \cdot e^2
!! \cdot [ (q+G)_y \cdot \text{Zstar}_{a,yx} \cdot \text{exp}
!! (-i\cdot(q+G)\cdot\text{tau}_a)) ]
!! / [ (q+G)_y \cdot \epsilon_yz \cdot (q+G)_z ] $$
!! a: atom index, x, y: Cartesian direction index
!!
!
!! Since QE uses Rydberg units, we multiply the e^2 = 2.0 factor.
!! The units of q and G are 2pi/a, so we need to divide dvlong by tpiba.
!!

6
PHonon/PH/dwfc.f90
View File

@ -27,7 +27,7 @@ SUBROUTINE dwfc (npw_ , igk_ , ik_ , icart_ , func, dfunc)
COMPLEX(DP), INTENT(IN) :: func(npwx)
COMPLEX(DP), INTENT(OUT) :: dfunc(npwx)
!
! Local variables
! ... local variables
!
INTEGER :: ig
REAL(DP) :: gvec, xk_aux
@ -55,7 +55,7 @@ SUBROUTINE d2wfc (npw_ , igk_ , ik_ , icart_ , jcart_, func, d2func)
!! This routine calculates the second derivative of a wave function
!! with respect to the position operator \(r_\text{icart}\) and
!! \(r_\text{jcart}\):
!! $$ d^2 \psi/d\text{icart}d\text{jcart} = \sum_G [-(k+G)_\text{icart}
!! $$ d^2 \psi/d\text{icart}\ d\text{jcart} = \sum_G [-(k+G)_\text{icart}
!! \cdot (-(k+G)_\text{gcart})] \psi(G) $$
!
USE kinds, ONLY : DP
@ -72,7 +72,7 @@ SUBROUTINE d2wfc (npw_ , igk_ , ik_ , icart_ , jcart_, func, d2func)
COMPLEX(DP), INTENT(IN) :: func(npwx)
COMPLEX(DP), INTENT(OUT) :: d2func(npwx)
!
! Local variables
! ... local variables
!
INTEGER :: ig
REAL(DP) :: gvec, xk_aux, gvec2, xk_aux2

18
PHonon/PH/dynmat_hub_bare.f90
View File

@ -9,19 +9,23 @@
SUBROUTINE dynmat_hub_bare
!---------------------------------------------------------------------
!! DFPT+U: This routine does two tasks:
!! 1) Computes d2ns_bare (very slow) or reads it;
!! 1) Computes \(\text{d2ns_bare}\) (very slow) or reads it;
!! 2) Adds to the dynamical matrix the bare part due the Hubbard term:
!! $$ \sum_\text{nah} \text{Hubbard_U}(\text{nah}) \sum_{is, m1, m2} [
!! (0.5\delta_{m1m2} - \text{ns}(m1,m2,is,\text{nah})) \cdot \text{d2ns_bare}
!! - \text{conjg}(\text{dnsbare}( m1, m2, is, \text{nah}, \text{icart},na))
!! \cdot \text{dnsbare}( m1, m2, is, nah, jcart,nap) ] $$
!! \begin{equation}\notag
!! \begin{split}
!! \sum_\text{nah} \text{Hubbard_U}&(\text{nah}) \sum_{is, m1, m2} [
!! (0.5\delta_{m1m2} - \text{ns}(m1,m2,is,\text{nah})) \cdot \text{d2ns_bare}\\
!! &- \text{conjg}(\text{dnsbare}( m1, m2, is, \text{nah}, \text{icart},\text{na}))
!! \cdot \text{dnsbare}( m1, m2, is, \text{nah}, \text{jcart},\text{nap}) ]
!! \end{split}
!! \end{equation}
!
!! Addition of the J0 term:
!! $$ \sum_\text{nah} \text{Hubbard_J0}(\text{nah}) \sum_{is, m1, m2} [
!! \text{ns}(m1,m2,-is,\text{nah}) \cdot \text{d2ns_bare}(m1, m2, is,
!! \text{nah}, \text{nap_jcart}, \text{na_icart})
!! + \text{conjg}(\text{dnsbare}( m1, m2, is, \text{nah}, \text{icart},na))\cdot
!! \text{dnsbare}( m1, m2, -is, nah, jcart,nap) ] $$
!! + \text{conjg}(\text{dnsbare}(m1,m2,is,\text{nah}, \text{icart},\text{na}))\cdot
!! \text{dnsbare}(m1,m2,-is, \text{nah}, \text{jcart},\text{nap}) ] $$
!
!! Written by A. Floris.
!! Modified by I. Timrov (01.10.2018).

7
PHonon/PH/dynmat_hub_scf.f90
View File

@ -13,9 +13,10 @@ SUBROUTINE dynmat_hub_scf (irr, nu_i0, nper)
!! It adds the scf derivative of the Hubbard energy (terms 1 and 2). Moreover,
!! it adds 3 terms due to the orthogonality constraints in the USPP formalism (terms 3,4,5).
!! Note, the orthogonality terms 4 and 5 DO NOT follow simply from the 2nd derivative of
!! the Hubbard energy \(\text{E_U}\); they stem from the coupling of the \(-\epsilon\dS\d\lambda\) in
!! the USPP forces with the part of d\psi\d\mu projected onto the occupied manifold
!! (the variable dpsi stems from the scf linear system and lives in the unoccupied manifold).
!! the Hubbard energy \(\text{E_U}\); they stem from the coupling of the \(-\epsilon\ dS\ d\lambda\) in
!! the USPP forces with the part of \(d\psi\ d\mu\) projected onto the occupied manifold
!! (the variable dpsi stems from the scf linear system and lives in the unoccupied manifold).
!! See comments in the source code for other details on the implemented terms.
!
! Terms implemented:
! dyn_hub_scf (ipert, imode)

2
PHonon/PH/gmressolve_all.f90
View File

@ -11,7 +11,7 @@ subroutine gmressolve_all (h_psi, cg_psi, e, d0psi, dpsi, h_diag, &
ndmx, ndim, ethr, ik, kter, conv_root, anorm, nbnd, m)
!----------------------------------------------------------------------
!! Iterative solution of the linear system by GMRES(m) method:
!! $$ (h - e + Q) \cdot \(\text{dpsi} = \(\text{d0psi}\)\ \ (1) $$
!! $$ (h - e + Q) \cdot \text{dpsi} = \text{d0psi}\ \ (1) $$
!! where \(h\) is a complex hermitian matrix, \(e\) is a
!! complex scalar and \(\text{dpsi}\) and \(\text{d0psi}\) are
!! complex vectors.

4
PHonon/PH/init_representations.f90
View File

@ -113,8 +113,8 @@ subroutine initialize_grid_variables()
!----------------------------------------------------------------------
!! This subroutine initializes the grid variables by reading the
!! modes from file. It uses the routine check_if_partial_dyn to
!! set the modes to compute according to \(\text{start_irr},
!! \(\text{last_irr} or \(\text{modenum}\) and \(\text{ifat}\) flags.
!! set the modes to compute according to \(\text{start_irr}\),
!! \(\text{last_irr}\) or \(\text{modenum}\) and \(\text{ifat}\) flags.
!
USE kinds, ONLY : DP
USE ions_base, ONLY : nat

18
PHonon/PH/phq_setup.f90
View File

@ -11,26 +11,26 @@ subroutine phq_setup
!! This subroutine prepares several variables which are needed in the
!! \(\texttt{phonon}\) program:
!! 1) computes the total local potential (external+scf) on the smooth
!! grid to be used in h_psi and similia;
!! grid to be used in \(\texttt{h_psi}\) and similia;
!! 2) computes the local magnetization (if necessary);
!! 3) computes dmuxc (with GC if needed);
!! 4) set the inverse of every matrix invs;
!! 5) for metals sets the occupied bands;
!! 6) computes alpha_pv;
!! 6) computes \(\text{alpha_pv}\);
!! 7) computes the variables needed to pass to the pattern representation:
!! \(\text{u}\): the patterns;
!! \(\text{t}\): the matrices of the small group of q on the pattern basis;
!! \(\text{tmq}\): the matrix of the symmetry which sends \(q\rightarrow-q+G}\);
!! \(\text{tmq}\): the matrix of the symmetry which sends \(q\rightarrow -q+G\);
!! \(\text{gi}\): the G associated to each symmetry operation;
!! \(\text{gimq}\): the G of the \(q\rigtharrow-q+G}\) symmetry;
!! \(\text{nsymq}\): the order of the small group of q;
!! \(\text{irotmq}\): the index of the \(q\rightarrow-q+G}\) symmetry;
!! \(\text{gimq}\): the G of the \(q\rightarrow -q+G\) symmetry;
!! \(\text{nsymq}\): the order of the small group of \(q\);
!! \(\text{irotmq}\): the index of the \(q\rightarrow -q+G\) symmetry;
!! \(\text{nirr}\): the number of irreducible representation;
!! \(\text{npert}\): the dimension of each irreducible representation;
!! \(\text{nmodes}\): the number of modes;
!! \(\text{minus_q}\): true if there is a symmetry sending \(q\rightarrow-q+G}\);
!! 8) for testing purposes it sets ubar;
!! 9) set the variables needed to deal with nlcc;
!! \(\text{minus_q}\): true if there is a symmetry sending \(q\rightarrow -q+G\);
!! 8) for testing purposes it sets \(\text{ubar}\);
!! 9) set the variables needed to deal with \(\text{nlcc}\);
!! 10) set the variables needed for the partial computation
! of the dynamical matrix.
!

2
PHonon/PH/q2qstar_ph.f90
View File

@ -11,7 +11,7 @@ subroutine q2qstar_ph (dyn, at, bg, nat, nsym, s, invs, irt, rtau, &
!-----------------------------------------------------------------------
!! Generates the dynamical matrices for the star of q and writes them on
!! disk for later use.
!! If there is a symmetry operation such that \(q \rigtharrow -q+G \) then
!! If there is a symmetry operation such that \(q \rightarrow -q+G \) then
!! imposes on dynamical matrix those conditions related to time reversal
!! symmetry.
!