mirror of https://gitlab.com/QEF/q-e.git
Temporary bug fix. A slower irrek is used in the noncollinear magnetic case,
because the old one seems to be wrong in some cases when some operations require time reversal. This solves bug number 68 indicated by Cezary Sliwa. git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@8231 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
parent
a8b2efa281
commit
cc23be5722
173
PW/irrek.f90
173
PW/irrek.f90
|
@ -17,6 +17,8 @@ subroutine irreducible_BZ (nrot, s, nsym, minus_q, at, bg, npk, nks, &
|
|||
! of the point group of the Bravais lattice.
|
||||
!
|
||||
USE kinds, only : DP
|
||||
USE noncollin_module, ONLY : noncolin
|
||||
USE spin_orb, ONLY : domag
|
||||
implicit none
|
||||
!
|
||||
integer, intent(in) :: nrot, nsym, npk, s(3,3,48), t_rev(48)
|
||||
|
@ -46,14 +48,20 @@ subroutine irreducible_BZ (nrot, s, nsym, minus_q, at, bg, npk, nks, &
|
|||
!
|
||||
! Find the coset in the point group of the Bravais lattice
|
||||
!
|
||||
sym(1:nsym) = .true.
|
||||
sym(nsym+1:)= .false.
|
||||
call coset (nrot, table, sym, nsym, irg)
|
||||
!
|
||||
! here we set the k-points in the irreducible wedge of the point grou
|
||||
! of the crystal
|
||||
!
|
||||
call irrek (at, bg, nrot, invs, nsym, irg, minus_q, npk, nks, xk, wk, t_rev)
|
||||
IF (noncolin.AND.domag) THEN
|
||||
call irrek_nc(at, bg, nrot, invs, nsym, irg, minus_q, npk, nks, xk, &
|
||||
wk, t_rev)
|
||||
ELSE
|
||||
sym(1:nsym) = .true.
|
||||
sym(nsym+1:)= .false.
|
||||
call coset (nrot, table, sym, nsym, irg)
|
||||
!
|
||||
! here we set the k-points in the irreducible wedge of the point grou
|
||||
! of the crystal
|
||||
!
|
||||
call irrek (at, bg, nrot, invs, nsym, irg, minus_q, npk, nks, xk, &
|
||||
wk, t_rev)
|
||||
ENDIF
|
||||
!
|
||||
return
|
||||
!
|
||||
|
@ -204,3 +212,152 @@ subroutine irrek (at, bg, nrot, invs, nsym, irg, minus_q, npk, &
|
|||
!
|
||||
return
|
||||
end subroutine irrek
|
||||
!-----------------------------------------------------------------------
|
||||
subroutine irrek_nc (at, bg, nrot, invs, nsym, irg, minus_q, npk, &
|
||||
nks, xk, wk, t_rev)
|
||||
!-----------------------------------------------------------------------
|
||||
!
|
||||
! Given a set of special points in the Irreducible Wedge of some
|
||||
! group, finds the equivalent special points in the IW of one of
|
||||
! its subgroups.
|
||||
!
|
||||
USE kinds, only : DP
|
||||
implicit none
|
||||
!
|
||||
integer, intent(in) :: npk, nrot, nsym, invs (3, 3, 48), irg (nrot)
|
||||
! maximum number of special points
|
||||
! order of the parent point group
|
||||
! order of the subgroup
|
||||
! inverse of the elements of the symmetry group
|
||||
! partition of the elements of the symmetry group into left cosets,
|
||||
! as given by SUBROUTINE COSET
|
||||
integer, intent(inout) :: nks
|
||||
! number of special points
|
||||
integer, intent(in) :: t_rev(48)
|
||||
real(DP), intent(in) :: at (3, 3), bg (3, 3)
|
||||
! basis vectors of the Bravais and reciprocal lattice
|
||||
real(DP), intent(inout) :: xk (3, npk), wk (npk)
|
||||
! special points and weights
|
||||
logical, intent(in) :: minus_q
|
||||
! .true. if symmetries q = -q+G are acceptable
|
||||
!
|
||||
! here the local variables
|
||||
!
|
||||
integer :: nks0, jk, kpol, irot, jrot, isym, ik, iks, start_k
|
||||
! nks0: used to save the initial number of k-points
|
||||
! ncos: total number of cosets
|
||||
real(DP) :: xkg (3), xks (3), xkn(3), one, xk_new(3,npk), wk_new(npk), &
|
||||
xk_cart(3)
|
||||
! coordinates of the k point in crystal axis
|
||||
! coordinates of the rotated k point
|
||||
! weight of each coset
|
||||
! buffer which contains the weight of k points
|
||||
! total weight of k-points
|
||||
logical :: satm
|
||||
! true if equivalent point found
|
||||
|
||||
nks0 = nks
|
||||
nks=0
|
||||
start_k=0
|
||||
DO jk = 1, nks0
|
||||
!
|
||||
! The k point is first computed in crystal axis
|
||||
!
|
||||
! xkg are the components of xk in the crystal base
|
||||
xkg (:) = at (1, :) * xk (1, jk) + &
|
||||
at (2, :) * xk (2, jk) + &
|
||||
at (3, :) * xk (3, jk)
|
||||
!
|
||||
! Then it is rotated with each symmetry of the global group.
|
||||
!
|
||||
DO irot = 1, nrot
|
||||
xks (:) = invs (:, 1, irot) * xkg (1) + &
|
||||
invs (:, 2, irot) * xkg (2) + &
|
||||
invs (:, 3, irot) * xkg (3)
|
||||
!
|
||||
! Now check if there is an operation of the subgroup that
|
||||
! makes xks equivalent to some other already found k point
|
||||
!
|
||||
DO jrot=1,nsym
|
||||
xkn (:) = invs (:, 1, jrot) * xks (1) + &
|
||||
invs (:, 2, jrot) * xks (2) + &
|
||||
invs (:, 3, jrot) * xks (3)
|
||||
IF (t_rev(jrot)==1) xkn =-xkn
|
||||
DO ik = start_k+1, nks
|
||||
satm = abs (xk_new (1, ik) - xkn (1) - &
|
||||
nint (xk_new (1, ik) - xkn (1) ) ) < 1.0d-5 .and. &
|
||||
abs (xk_new (2, ik) - xkn (2) - &
|
||||
nint (xk_new (2, ik) - xkn (2) ) ) < 1.0d-5 .and. &
|
||||
abs (xk_new (3, ik) - xkn (3) - &
|
||||
nint (xk_new (3, ik) - xkn (3) ) ) < 1.0d-5
|
||||
!
|
||||
! .... or equivalent to minus each other when minus_q=.t.
|
||||
!
|
||||
if (minus_q) satm = satm .or. &
|
||||
abs (xk_new (1, ik) + xkn (1) - &
|
||||
nint (xk_new (1, ik) + xkn (1) ) ) < 1.0d-5 .and. &
|
||||
abs (xk_new (2, ik) + xkn (2) - &
|
||||
nint (xk_new (2, ik) + xkn (2) ) ) < 1.0d-5 .and. &
|
||||
abs (xk_new (3, ik) + xkn (3) - &
|
||||
nint (xk_new (3, ik) + xkn (3) ) ) < 1.0d-5
|
||||
IF ( satm ) THEN
|
||||
wk_new(ik) = wk_new(ik) + wk(jk)
|
||||
GOTO 100
|
||||
ENDIF
|
||||
END DO
|
||||
END DO
|
||||
nks=nks+1
|
||||
IF (nks > npk) CALL errore('irrek_nc','too many k points',1)
|
||||
xk_new(:,nks)=xks
|
||||
wk_new(nks)=wk(jk)
|
||||
100 CONTINUE
|
||||
ENDDO
|
||||
start_k=nks
|
||||
ENDDO
|
||||
!
|
||||
! The order of the original k points is preserved
|
||||
!
|
||||
iks=nks0
|
||||
DO ik = 1, nks
|
||||
!
|
||||
! for each new k point found, check if it was in the original list
|
||||
!
|
||||
DO jk=1, nks0
|
||||
xkg (:) = at (1, :) * xk (1, jk) + &
|
||||
at (2, :) * xk (2, jk) + &
|
||||
at (3, :) * xk (3, jk)
|
||||
satm = abs (xk_new (1, ik) - xkg (1) - &
|
||||
nint (xk_new (1, ik) - xkg (1) ) ) < 1.0d-5 .and. &
|
||||
abs (xk_new (2, ik) - xkg (2) - &
|
||||
nint (xk_new (2, ik) - xkg (2) ) ) < 1.0d-5 .and. &
|
||||
abs (xk_new (3, ik) - xkg (3) - &
|
||||
nint (xk_new (3, ik) - xkg (3) ) ) < 1.0d-5
|
||||
IF (satm) THEN
|
||||
!
|
||||
! If it was, just update the weight
|
||||
!
|
||||
wk(jk)=wk_new(ik)
|
||||
goto 200
|
||||
ENDIF
|
||||
ENDDO
|
||||
!
|
||||
! If it was not, bring xk_new in cartesian coodinates and copy it in the
|
||||
! first free place available
|
||||
!
|
||||
iks=iks+1
|
||||
xk_cart (:) = bg (:, 1) * xk_new (1, ik) + &
|
||||
bg (:, 2) * xk_new (2, ik) + &
|
||||
bg (:, 3) * xk_new (3, ik)
|
||||
xk(:,iks)=xk_cart(:)
|
||||
wk(iks)=wk_new(ik)
|
||||
200 CONTINUE
|
||||
ENDDO
|
||||
IF (iks /= nks ) CALL errore('irrek_nc','Internal problem with k points',1)
|
||||
!
|
||||
! normalize weights to one
|
||||
!
|
||||
one = SUM (wk(1:nks))
|
||||
IF ( one > 0.d0 ) wk(1:nks) = wk(1:nks) / one
|
||||
!
|
||||
RETURN
|
||||
END SUBROUTINE irrek_nc
|
||||
|
|
Loading…
Reference in New Issue