TABLE OF CONTENTS


ABINIT/cgtk_rotate [ Functions ]

[ Top ] [ Functions ]

NAME

  cgtk_rotate

FUNCTION

INPUTS

  isym
  itimrev
  nspinor
  ndat
  npw1,npw2
  istwf1,istwf2
  cryst<crystal_t>
  shiftg(3)
  kg1(3,npw1), kg2(3,npw2)
  work_ngfft(18)
  kpt1(3)
  cg1(2,npw1,nspinor,ndat)

OUTPUT

  cg2(2,npw2,nspinor,ndat)
  work(2,work_ngfft(4),work_ngfft(5),work_ngfft(6)) !*ndat),

PARENTS

      m_phgamma,m_sigmaph,m_wfk

NOTES

  Inspired to wfconv.

CHILDREN

      getph,getspinrot,mati3inv,ph1d3d,sphere

SOURCE

 85 subroutine cgtk_rotate(cryst,kpt1,isym,itimrev,shiftg,nspinor,ndat,&
 86   npw1,kg1,npw2,kg2,istwf1,istwf2,cg1,cg2,work_ngfft,work)
 87 
 88 
 89 !This section has been created automatically by the script Abilint (TD).
 90 !Do not modify the following lines by hand.
 91 #undef ABI_FUNC
 92 #define ABI_FUNC 'cgtk_rotate'
 93 !End of the abilint section
 94 
 95  implicit none
 96 
 97 !Arguments ------------------------------------
 98 !scalars
 99  integer,intent(in) :: isym,itimrev,nspinor,ndat,npw1,npw2,istwf1,istwf2
100  type(crystal_t),intent(in) :: cryst
101 !arrays
102  integer,intent(in) :: shiftg(3),kg1(3,npw1),kg2(3,npw2)
103  integer,intent(in) :: work_ngfft(18)
104  real(dp),intent(in) :: kpt1(3)
105  real(dp),intent(in) :: cg1(2,npw1,nspinor,ndat)
106  real(dp),intent(out) :: cg2(2,npw2,nspinor,ndat)
107  real(dp),intent(out) :: work(2,work_ngfft(4),work_ngfft(5),work_ngfft(6)) !*ndat) for threads?
108 
109 !Local variables ------------------------------
110 !scalars
111  integer,parameter :: tobox=1,tosph=-1,me_g0=1
112  integer :: n1,n2,n3,n4,n5,n6,ipw,idat,isp
113  real(dp) :: arg,ar,ai,bi,br,spinrots,spinrotx,spinroty,spinrotz
114  logical :: have_phase
115 !arrays
116  integer,parameter :: no_shift(3)=0,atindx(1)=1
117  integer :: symm(3,3),symrel_conv(3,3)
118  real(dp) :: phktnons(2,1),tnons_conv(3),spinrot(4)
119  real(dp),allocatable :: phase1d(:,:),phase3d(:,:),wavef1(:,:)
120 
121 !************************************************************************
122 
123  n1=work_ngfft(1); n2=work_ngfft(2); n3=work_ngfft(3)
124  n4=work_ngfft(4); n5=work_ngfft(5); n6=work_ngfft(6)
125 
126  symrel_conv = cryst%symrel(:,:,isym)
127  call mati3inv(symrel_conv, symm)
128  tnons_conv = cryst%tnons(:,isym) !
129  have_phase = sum(tnons_conv**2) > tol8
130 
131  ! Compute rotation in spinor space
132  if (nspinor == 2) call getspinrot(cryst%rprimd,spinrot,symrel_conv)
133  if (itimrev == 1) symm = -symm
134 
135  ! Need to compute phase factors associated with nonsymmorphic translations?
136  if (have_phase) then
137    ABI_MALLOC(phase3d,(2,npw1))
138    ABI_MALLOC(phase1d,(2,(2*n1+1)+(2*n2+1)+(2*n3+1)))
139    ! Although the routine getph is originally written for atomic phase factors, it does precisely what we want
140    call getph(atindx,1,n1,n2,n3,phase1d,tnons_conv)
141 
142    arg = two_pi*(kpt1(1)*tnons_conv(1)+ kpt1(2)*tnons_conv(2)+ kpt1(3)*tnons_conv(3))
143    phktnons(1,1)=cos(arg)
144    phktnons(2,1)=sin(arg)
145    ! Convert 1D phase factors to 3D phase factors exp(i 2 pi (k+G).tnons )
146    call ph1d3d(1,1,kg1,1,1,npw1,n1,n2,n3,phktnons,phase1d,phase3d)
147  end if
148 
149  ABI_MALLOC(wavef1, (2,npw1))
150 
151  do idat=1,ndat
152 
153    do isp=1,nspinor
154      wavef1 = cg1(:,:,isp,idat)
155 
156      if (have_phase) then
157        ! Multiply by phase factors due to nonsymmorphic translations.
158        do ipw=1,npw1
159          ar = phase3d(1,ipw)*wavef1(1,ipw)-phase3d(2,ipw)*wavef1(2,ipw)
160          ai = phase3d(2,ipw)*wavef1(1,ipw)+phase3d(1,ipw)*wavef1(2,ipw)
161          wavef1(1,ipw)=ar
162          wavef1(2,ipw)=ai
163        end do
164      end if
165 
166      ! Take into account time-reversal symmetry, if needed, in the scalar case
167      if (itimrev == 1 .and. nspinor == 1) wavef1(2,:npw1) = -wavef1(2,:npw1)
168 
169      ! Insert wavef1 in work array.
170      call sphere(wavef1,1,npw1,work,n1,n2,n3,n4,n5,n6,kg1,istwf1,tobox,me_g0,no_shift,identity_3d,one)
171 
172      ! Apply rotation and extract data.
173      call sphere(cg2(:,:,isp,idat),1,npw2,work,n1,n2,n3,n4,n5,n6,kg2,istwf2,tosph,me_g0,shiftg,symm,one)
174    end do ! isp
175 
176    if (nspinor == 2) then
177      if (itimrev == 1) then
178        ! Take care of time-reversal symmetry, if needed
179        ! Exchange spin-up and spin-down. Make complex conjugate of one component,
180        ! and change sign of other component
181        do ipw=1,npw2
182          ! Here, change sign of real part
183          ar = -cg2(1,ipw,1,idat)
184          ai =  cg2(2,ipw,1,idat)
185          ! Here, change sign of imaginary part
186          cg2(1,ipw,1,idat) =  cg2(1,ipw,2,idat)
187          cg2(2,ipw,1,idat) = -cg2(2,ipw,2,idat)
188          cg2(1,ipw,2,idat) = ar
189          cg2(2,ipw,2,idat) = ai
190        end do
191      end if ! itimrev==1
192 
193      ! Rotation in spinor space (see also wfconv)
194      spinrots = spinrot(1)
195      spinrotx = spinrot(2)
196      spinroty = spinrot(3)
197      spinrotz = spinrot(4)
198      do ipw=1,npw2
199        ar = cg2(1,ipw,1,idat)
200        ai = cg2(2,ipw,1,idat)
201        br = cg2(1,ipw,2,idat)
202        bi = cg2(2,ipw,2,idat)
203        cg2(1,ipw,1,idat) =  spinrots*ar-spinrotz*ai +spinroty*br-spinrotx*bi
204        cg2(2,ipw,1,idat) =  spinrots*ai+spinrotz*ar +spinroty*bi+spinrotx*br
205        cg2(1,ipw,2,idat) = -spinroty*ar-spinrotx*ai +spinrots*br+spinrotz*bi
206        cg2(2,ipw,2,idat) = -spinroty*ai+spinrotx*ar +spinrots*bi-spinrotz*br
207      end do
208    end if
209 
210  end do ! idat
211 
212  ABI_FREE(wavef1)
213 
214  if (have_phase) then
215    ABI_FREE(phase3d)
216    ABI_FREE(phase1d)
217  end if
218 
219 end subroutine cgtk_rotate

ABINIT/m_cgtk [ Modules ]

[ Top ] [ Modules ]

NAME

  m_cgtk

FUNCTION

COPYRIGHT

  Copyright (C) 2008-2018 ABINIT group (MG)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_cgtk
28 
29  use defs_basis
30  use m_errors
31  use m_abicore
32 
33  use m_symtk,     only : mati3inv
34  use m_geometry,  only : getspinrot
35  use m_crystal,   only : crystal_t
36  use m_fftcore,   only : sphere
37  use m_kg,        only : ph1d3d, getph
38 
39  implicit none
40 
41  private