TABLE OF CONTENTS


ABINIT/m_cplxtools [ Modules ]

[ Top ] [ Modules ]

NAME

  m_cplxtools

FUNCTION

 This module defines helper functions to operate on complex arrays (mainly used in the GW code)

COPYRIGHT

 Copyright (C) 1992-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 .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

NOTES

 1) The convention about names of interfaced routine is: cplx_<name>,
    where <name> is equal to the name of the standard BLAS routine


m_cplxtools/cplx_addtorho_dpc [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_addtorho_dpc

FUNCTION

  Add |ur|**2 to the ground-states density rho.
    rho = rho + weight_r * |ur|**2

INPUTS

  nx,ny,nz=physical dimension of the FFT box.
  ldx,ldy,ldz=leading dimensions of the arrays.
  ndat=number of contributions to accumulate.
  weight_r=weight used for the accumulation of the density in real space
  ur(ldx,ldy,ldz*ndat)=wavefunctions in real space

SIDE EFFECTS

  rho(ldx,ldy,ldz) = contains the input density at input,
                     modified in input with the contribution gived by ur.

PARENTS

CHILDREN

SOURCE

1209 subroutine cplx_addtorho_dpc(nx,ny,nz,ldx,ldy,ldz,ndat,weight_r,ur,rho)
1210 
1211 
1212 !This section has been created automatically by the script Abilint (TD).
1213 !Do not modify the following lines by hand.
1214 #undef ABI_FUNC
1215 #define ABI_FUNC 'cplx_addtorho_dpc'
1216 !End of the abilint section
1217 
1218  implicit none
1219 
1220 !Arguments ------------------------------------
1221 !scalars
1222  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat
1223  real(dp),intent(in) :: weight_r
1224 !arrays
1225  complex(dpc),intent(in) :: ur(ldx*ldy*ldz*ndat)
1226  real(dp),intent(inout) :: rho(ldx*ldy*ldz)
1227 
1228 !Local variables-------------------------------
1229 !scalars
1230  integer :: ix,iy,iz,dat,ifft,ldxyz,pad_box,padz,pady
1231 
1232 ! *************************************************************************
1233 
1234  ldxyz = ldx*ldy*ldz
1235 
1236  if (ndat==1) then
1237 !$OMP PARALLEL DO PRIVATE(padz, pady, ifft)
1238    do iz=1,nz
1239      padz = (iz-1)*ldx*ldy
1240      do iy=1,ny
1241        pady = (iy-1)*ldx
1242        do ix=1,nx
1243          ifft = ix + pady + padz
1244          rho(ifft) = rho(ifft) + weight_r * (REAL(ur(ifft))**2 + AIMAG(ur(ifft))**2)
1245        end do
1246      end do
1247    end do
1248 
1249  else
1250 ! It would be nice to use $OMP PARALLEL DO REDUCTION(+:rho)
1251 ! but it's risky as the private rho is allocated on the stack of the thread.
1252 !$OMP PARALLEL PRIVATE(pad_box, padz, pady, ifft)
1253    do dat=1,ndat
1254      pad_box = (dat-1)*ldxyz
1255 !$OMP DO
1256      do iz=1,nz
1257        padz = (iz-1)*ldx*ldy
1258        do iy=1,ny
1259          pady = (iy-1)*ldx
1260          do ix=1,nx
1261            ifft = ix + pady + padz
1262            rho(ifft) = rho(ifft) + weight_r * (REAL(ur(ifft+pad_box)**2 + AIMAG(ur(ifft+pad_box))**2))
1263          end do
1264        end do
1265      end do
1266 !$OMP END DO NOWAIT
1267    end do
1268 !$OMP END PARALLEL
1269  end if
1270 
1271 end subroutine cplx_addtorho_dpc

m_cplxtools/cplx_box2gsph_dpc [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_box2gsph_dpc

FUNCTION

   Transfer data from the FFT box to the G-sphere. Target DPC complex array.

INPUTS

  nx,ny,nz=physical dimension of the FFT box.
  ldx,ldy,ldz=Logical dimensions of the arrays.
  ndat=number of data in iarrbox
  npw_k=Number of planewaves in the G-sphere.
  kg_k(3,npw_k)=Reduced coordinates of the G-vectoes.
  iarrbox(ldx*ldy*ldz*ndat)=Complex Input arrays on the FFT box.
  [rscal] = Scaling factor

OUTPUT

  oarrsph(npw_k*ndat)=Complex Data defined on the G-sphere.

PARENTS

CHILDREN

SOURCE

602 subroutine cplx_box2gsph_dpc(nx,ny,nz,ldx,ldy,ldz,ndat,npw_k,kg_k,iarrbox,oarrsph,rscal)
603 
604 
605 !This section has been created automatically by the script Abilint (TD).
606 !Do not modify the following lines by hand.
607 #undef ABI_FUNC
608 #define ABI_FUNC 'cplx_box2gsph_dpc'
609 !End of the abilint section
610 
611  implicit none
612 
613 !Arguments ------------------------------------
614 !scalars
615  integer,intent(in) :: npw_k,nx,ny,nz,ldx,ldy,ldz,ndat
616  real(dp),optional,intent(in) :: rscal
617 !arrays
618  integer,intent(in) :: kg_k(3,npw_k)
619  complex(dpc),intent(in) :: iarrbox(ldx*ldy*ldz*ndat)
620  complex(dpc),intent(out) :: oarrsph(npw_k*ndat)
621 
622 !Local variables-------------------------------
623 !scalars
624  integer :: ig,ix,iy,iz,dat,pad_sph,pad_box,ifft,ldxyz
625 
626 ! *************************************************************************
627 
628  ldxyz = ldx*ldy*ldz
629  if (.not. PRESENT(rscal)) then
630    !
631    if (ndat==1) then
632 !$OMP PARALLEL DO PRIVATE(ix,iy,iz,ifft)
633      do ig=1,npw_k
634        ix=kg_k(1,ig); if (ix<0) ix=ix+nx; ix=ix+1
635        iy=kg_k(2,ig); if (iy<0) iy=iy+ny; iy=iy+1
636        iz=kg_k(3,ig); if (iz<0) iz=iz+nz; iz=iz+1
637        ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
638        oarrsph(ig) = iarrbox(ifft)
639      end do
640    else
641 !$OMP PARALLEL DO PRIVATE(pad_sph,pad_box,ix,iy,iz,ifft)
642      do dat=1,ndat
643        pad_sph = (dat-1)*npw_k
644        pad_box = (dat-1)*ldxyz
645        do ig=1,npw_k
646          ix=kg_k(1,ig); if (ix<0) ix=ix+nx; ix=ix+1
647          iy=kg_k(2,ig); if (iy<0) iy=iy+ny; iy=iy+1
648          iz=kg_k(3,ig); if (iz<0) iz=iz+nz; iz=iz+1
649          ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
650          oarrsph(ig+pad_sph) = iarrbox(ifft+pad_box)
651        end do
652      end do
653    end if
654    !
655  else
656    if (ndat==1) then
657 !$OMP PARALLEL DO PRIVATE(ix,iy,iz,ifft)
658      do ig=1,npw_k
659        ix=kg_k(1,ig); if (ix<0) ix=ix+nx; ix=ix+1
660        iy=kg_k(2,ig); if (iy<0) iy=iy+ny; iy=iy+1
661        iz=kg_k(3,ig); if (iz<0) iz=iz+nz; iz=iz+1
662        ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
663        oarrsph(ig) = iarrbox(ifft) * rscal
664      end do
665    else
666 !$OMP PARALLEL DO PRIVATE(pad_sph,pad_box,ix,iy,iz,ifft)
667      do dat=1,ndat
668        pad_sph = (dat-1)*npw_k
669        pad_box = (dat-1)*ldxyz
670        do ig=1,npw_k
671          ix=kg_k(1,ig); if (ix<0) ix=ix+nx; ix=ix+1
672          iy=kg_k(2,ig); if (iy<0) iy=iy+ny; iy=iy+1
673          iz=kg_k(3,ig); if (iz<0) iz=iz+nz; iz=iz+1
674          ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
675          oarrsph(ig+pad_sph) = iarrbox(ifft+pad_box) * rscal
676        end do
677      end do
678    end if
679  end if
680 
681 end subroutine cplx_box2gsph_dpc

m_cplxtools/cplx_box2gsph_spc [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_box2gsph_spc

FUNCTION

   Transfer data from the FFT box to the G-sphere. Target SPC complex array.

INPUTS

  nx,ny,nz=physical dimension of the FFT box.
  ldx,ldy,ldz=Logical dimensions of the arrays.
  ndat=number of data in iarrbox
  npw_k=Number of planewaves in the G-sphere.
  kg_k(3,npw_k)=Reduced coordinates of the G-vectoes.
  iarrbox(ldx*ldy*ldz*ndat)=Complex Input arrays on the FFT box.
  [rscal] = Scaling factor

OUTPUT

  oarrsph(npw_k*ndat)=Complex Data defined on the G-sphere.

PARENTS

CHILDREN

SOURCE

491 subroutine cplx_box2gsph_spc(nx,ny,nz,ldx,ldy,ldz,ndat,npw_k,kg_k,iarrbox,oarrsph,rscal)
492 
493 
494 !This section has been created automatically by the script Abilint (TD).
495 !Do not modify the following lines by hand.
496 #undef ABI_FUNC
497 #define ABI_FUNC 'cplx_box2gsph_spc'
498 !End of the abilint section
499 
500  implicit none
501 
502 !Arguments ------------------------------------
503 !scalars
504  integer,intent(in) :: npw_k,nx,ny,nz,ldx,ldy,ldz,ndat
505  real(sp),optional,intent(in) :: rscal
506 !arrays
507  integer,intent(in) :: kg_k(3,npw_k)
508  complex(spc),intent(in) :: iarrbox(ldx*ldy*ldz*ndat)
509  complex(spc),intent(out) :: oarrsph(npw_k*ndat)
510 
511 !Local variables-------------------------------
512 !scalars
513  integer :: ig,ix,iy,iz,dat,pad_sph,pad_box,ifft,ldxyz
514 
515 ! *************************************************************************
516 
517  ldxyz = ldx*ldy*ldz
518  if (.not. PRESENT(rscal)) then
519    !
520    if (ndat==1) then
521 !$OMP PARALLEL DO PRIVATE(ix,iy,iz,ifft)
522      do ig=1,npw_k
523        ix=kg_k(1,ig); if (ix<0) ix=ix+nx; ix=ix+1
524        iy=kg_k(2,ig); if (iy<0) iy=iy+ny; iy=iy+1
525        iz=kg_k(3,ig); if (iz<0) iz=iz+nz; iz=iz+1
526        ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
527        oarrsph(ig) = iarrbox(ifft)
528      end do
529    else
530 !$OMP PARALLEL DO PRIVATE(pad_sph,pad_box,ix,iy,iz,ifft)
531      do dat=1,ndat
532        pad_sph = (dat-1)*npw_k
533        pad_box = (dat-1)*ldxyz
534        do ig=1,npw_k
535          ix=kg_k(1,ig); if (ix<0) ix=ix+nx; ix=ix+1
536          iy=kg_k(2,ig); if (iy<0) iy=iy+ny; iy=iy+1
537          iz=kg_k(3,ig); if (iz<0) iz=iz+nz; iz=iz+1
538          ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
539          oarrsph(ig+pad_sph) = iarrbox(ifft+pad_box)
540        end do
541      end do
542    end if
543    !
544  else
545    if (ndat==1) then
546 !$OMP PARALLEL DO PRIVATE(ix,iy,iz,ifft)
547      do ig=1,npw_k
548        ix=kg_k(1,ig); if (ix<0) ix=ix+nx; ix=ix+1
549        iy=kg_k(2,ig); if (iy<0) iy=iy+ny; iy=iy+1
550        iz=kg_k(3,ig); if (iz<0) iz=iz+nz; iz=iz+1
551        ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
552        oarrsph(ig) = iarrbox(ifft) * rscal
553      end do
554    else
555 !$OMP PARALLEL DO PRIVATE(pad_sph,pad_box,ix,iy,iz,ifft)
556      do dat=1,ndat
557        pad_sph = (dat-1)*npw_k
558        pad_box = (dat-1)*ldxyz
559        do ig=1,npw_k
560          ix=kg_k(1,ig); if (ix<0) ix=ix+nx; ix=ix+1
561          iy=kg_k(2,ig); if (iy<0) iy=iy+ny; iy=iy+1
562          iz=kg_k(3,ig); if (iz<0) iz=iz+nz; iz=iz+1
563          ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
564          oarrsph(ig+pad_sph) = iarrbox(ifft+pad_box) * rscal
565        end do
566      end do
567    end if
568  end if
569 
570 end subroutine cplx_box2gsph_spc

m_cplxtools/cplx_filter [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_filter

FUNCTION

  Set all the elements of x to zero where mask is .TRUE.

INPUTS

  n=Specifies the number of elements in vectors x and y.
  mask(n)=Logical array.

SIDE EFFECTS

  x(n)=See description.

PARENTS

CHILDREN

SOURCE

171 subroutine cplx_filter(n, x, mask)
172 
173 
174 !This section has been created automatically by the script Abilint (TD).
175 !Do not modify the following lines by hand.
176 #undef ABI_FUNC
177 #define ABI_FUNC 'cplx_filter'
178 !End of the abilint section
179 
180  implicit none
181 
182 !Arguments ------------------------------------
183 !scalars
184  integer,intent(in) :: n
185 !arrays
186  complex(dpc),intent(inout) :: x(n)
187  logical,intent(in) :: mask(n)
188 
189 ! *************************************************************************
190 
191  where (mask)
192    x = czero
193  end where
194 
195 end subroutine cplx_filter

m_cplxtools/cplx_fromreal [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_fromreal

FUNCTION

  Convert a real array with (real,imag) part to a complex array

INPUTS

  n = Specifies the number of elements in ocplx
  ireal(2*n)=Input real array.

OUTPUT

  ocplx(n)=Output complex array

PARENTS

CHILDREN

SOURCE

117 subroutine cplx_fromreal(n,ireal,ocplx)
118 
119 
120 !This section has been created automatically by the script Abilint (TD).
121 !Do not modify the following lines by hand.
122 #undef ABI_FUNC
123 #define ABI_FUNC 'cplx_fromreal'
124 !End of the abilint section
125 
126  implicit none
127 
128 !Arguments ------------------------------------
129 !scalars
130  integer,intent(in) :: n
131 !arrays
132  real(dp),intent(in) :: ireal(2,n)
133  complex(dpc),intent(out) :: ocplx(n)
134 
135 !Local variables ------------------------------
136 !scalars
137  integer :: ii
138 
139 ! *************************************************************************
140 
141 !$OMP PARALLEL DO PRIVATE(ii)
142  do ii=1,n
143    ocplx(ii) = DCMPLX(ireal(1,ii),ireal(2,ii))
144  end do
145 
146 end subroutine cplx_fromreal

m_cplxtools/cplx_gsph2box_dpc [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

 cplx_gsph2box_dpc

FUNCTION

 Array iarrsph is defined in sphere with npw points. Insert iarrsph inside box
 of nx*ny*nz points to define array oarrbox for fft box. rest of oarrbox is filled with 0 s.
 targer: DPC complex arrays

INPUTS

 iarrsph(2,npw*ndat)= contains values for npw G vectors in basis sphere
 ndat=number of FFT to perform.
 npw=number of G vectors in basis at this k point
 oarrbox(2,ldx*ldy*ldz*ndat) = fft box
 nx,ny,nz=physical dimension of the box (oarrbox)
 ldx,ldy,ldz=memory dimension of oarrbox
 kg_k(3,npw)=integer coordinates of G vectors in basis sphere
 istwf_k=option parameter that describes the storage of wfs

OUTPUT

   oarrbox(ldx*ldy*ldz*ndat)

NOTES

 If istwf_k differs from 1, then special storage modes must be taken
 into account, for symmetric wavefunctions coming from k=(0 0 0) or other
 special k points.

PARENTS

CHILDREN

SOURCE

 896 subroutine cplx_gsph2box_dpc(nx,ny,nz,ldx,ldy,ldz,ndat,npw,istwf_k,kg_k,iarrsph,oarrbox)
 897 
 898 
 899 !This section has been created automatically by the script Abilint (TD).
 900 !Do not modify the following lines by hand.
 901 #undef ABI_FUNC
 902 #define ABI_FUNC 'cplx_gsph2box_dpc'
 903 !End of the abilint section
 904 
 905  implicit none
 906 
 907 !Arguments ------------------------------------
 908 !scalars
 909  integer,intent(in) :: istwf_k,nx,ny,nz,ldx,ldy,ldz,ndat,npw
 910 !arrays
 911  integer,intent(in) :: kg_k(3,npw)
 912  complex(dpc),intent(in) :: iarrsph(npw*ndat)
 913  complex(dpc),intent(out) :: oarrbox(ldx*ldy*ldz*ndat)
 914 
 915 !Local variables-------------------------------
 916 !scalars
 917  integer,parameter :: me_g0=1
 918  integer :: ix,ixinv,iy,iyinv,iz,izinv,dat,ipw,npwmin,pad_box,pad_sph,ifft,ifft_inv,ldxyz
 919  !character(len=500) :: msg
 920 !arrays
 921  integer,allocatable :: ixinver(:),iyinver(:),izinver(:)
 922 
 923 ! *************************************************************************
 924 
 925 !In the case of special k-points, invariant under time-reversal,
 926 !but not Gamma, initialize the inverse coordinates
 927 !Remember indeed that
 928 !u_k(G) = u_{k+G0}(G-G0); u_{-k}(G) = u_k(G)^*
 929 !and therefore:
 930 !u_{G0/2}(G) = u_{G0/2}(-G-G0)^*.
 931  if (istwf_k>=2) then
 932    ABI_MALLOC(ixinver,(nx))
 933    ABI_MALLOC(iyinver,(ny))
 934    ABI_MALLOC(izinver,(nz))
 935    if ( ANY(istwf_k==(/2,4,6,8/)) ) then
 936      ixinver(1)=1
 937      do ix=2,nx
 938        ixinver(ix)=nx+2-ix
 939      end do
 940    else
 941      do ix=1,nx
 942        ixinver(ix)=nx+1-ix
 943      end do
 944    end if
 945    if (istwf_k>=2 .and. istwf_k<=5) then
 946      iyinver(1)=1
 947      do iy=2,ny
 948        iyinver(iy)=ny+2-iy
 949      end do
 950    else
 951      do iy=1,ny
 952        iyinver(iy)=ny+1-iy
 953      end do
 954    end if
 955    if ( ANY(istwf_k==(/2,3,6,7/)) ) then
 956      izinver(1)=1
 957      do iz=2,nz
 958        izinver(iz)=nz+2-iz
 959      end do
 960    else
 961      do iz=1,nz
 962        izinver(iz)=nz+1-iz
 963      end do
 964    end if
 965  end if
 966 
 967  ldxyz = ldx*ldy*ldz
 968 
 969  if (istwf_k==1) then
 970 
 971 !$OMP PARALLEL DO PRIVATE(pad_sph,pad_box,ix,iy,iz,ifft)
 972    do dat=1,ndat
 973      pad_sph = (dat-1)*npw
 974      pad_box = (dat-1)*ldxyz
 975      oarrbox(1+pad_box:ldxyz+pad_box) = czero_dpc ! zero the sub-array
 976      do ipw=1,npw
 977        ix=kg_k(1,ipw); if (ix<0) ix=ix+nx; ix=ix+1
 978        iy=kg_k(2,ipw); if (iy<0) iy=iy+ny; iy=iy+1
 979        iz=kg_k(3,ipw); if (iz<0) iz=iz+nz; iz=iz+1
 980        ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
 981 #if defined __INTEL_COMPILER && defined HAVE_OPENMP
 982        if (ifft==0) then
 983          MSG_ERROR("prevent ifort+OMP from miscompiling this section on cronos")
 984        end if
 985 #endif
 986        oarrbox(ifft+pad_box) = iarrsph(ipw+pad_sph)
 987      end do
 988    end do
 989 
 990  else if (istwf_k>=2) then
 991    !
 992    npwmin=1
 993    if(istwf_k==2 .and. me_g0==1) then ! If gamma point, then oarrbox must be completed
 994      do dat=1,ndat
 995        pad_sph = (dat-1)*npw
 996        pad_box = (dat-1)*ldxyz
 997        oarrbox(1+pad_box) = REAL(iarrsph(1+pad_sph))
 998      end do
 999      npwmin=2
1000    end if
1001 
1002 !$OMP PARALLEL DO PRIVATE(pad_sph,pad_box,ix,ixinv,iy,iyinv,iz,izinv,ifft)
1003    do dat=1,ndat
1004      pad_sph = (dat-1)*npw
1005      pad_box = (dat-1)*ldxyz
1006      oarrbox(npwmin+pad_box:ldxyz+pad_box) = czero_dpc
1007      do ipw=npwmin,npw
1008        ix=kg_k(1,ipw); if(ix<0)ix=ix+nx; ix=ix+1
1009        iy=kg_k(2,ipw); if(iy<0)iy=iy+ny; iy=iy+1
1010        iz=kg_k(3,ipw); if(iz<0)iz=iz+nz; iz=iz+1
1011        ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
1012        ! Construct the coordinates of -k-G
1013        ixinv=ixinver(ix); iyinv=iyinver(iy); izinv=izinver(iz)
1014        ifft_inv = ixinv + (iyinv-1)*ldx + (izinv-1)*ldx*ldy
1015 #if defined __INTEL_COMPILER && defined HAVE_OPENMP
1016        if (ifft==0 .or. ifft_inv==0) then
1017          MSG_ERROR("prevent ifort+OMP from miscompiling this section on cronos")
1018        end if
1019 #endif
1020        oarrbox(ifft    +pad_box) =        iarrsph(ipw+pad_sph)
1021        oarrbox(ifft_inv+pad_box) = DCONJG(iarrsph(ipw+pad_sph))
1022      end do
1023    end do
1024    !
1025  else
1026    MSG_ERROR("Wrong istwfk")
1027  end if
1028 
1029  if (istwf_k>=2) then
1030    ABI_FREE(ixinver)
1031    ABI_FREE(iyinver)
1032    ABI_FREE(izinver)
1033  end if
1034 
1035 end subroutine cplx_gsph2box_dpc

m_cplxtools/cplx_gsph2box_spc [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

 cplx_gsph2box_spc

FUNCTION

 Array iarrsph is defined in sphere with npw points. Insert iarrsph inside box
 of nx*ny*nz points to define array oarrbox for fft box. rest of oarrbox is filled with 0 s.
 targer: SPC complex arrays

INPUTS

 iarrsph(2,npw*ndat)= contains values for npw G vectors in basis sphere
 ndat=number of FFT to perform.
 npw=number of G vectors in basis at this k point
 oarrbox(2,ldx*ldy*ldz*ndat) = fft box
 nx,ny,nz=physical dimension of the box (oarrbox)
 ldx,ldy,ldz=memory dimension of oarrbox
 kg_k(3,npw)=integer coordinates of G vectors in basis sphere
 istwf_k=option parameter that describes the storage of wfs

OUTPUT

   oarrbox(ldx*ldy*ldz*ndat)

NOTES

 If istwf_k differs from 1, then special storage modes must be taken
 into account, for symmetric wavefunctions coming from k=(0 0 0) or other
 special k points.

PARENTS

CHILDREN

SOURCE

719 subroutine cplx_gsph2box_spc(nx,ny,nz,ldx,ldy,ldz,ndat,npw,istwf_k,kg_k,iarrsph,oarrbox)
720 
721 
722 !This section has been created automatically by the script Abilint (TD).
723 !Do not modify the following lines by hand.
724 #undef ABI_FUNC
725 #define ABI_FUNC 'cplx_gsph2box_spc'
726 !End of the abilint section
727 
728  implicit none
729 
730 !Arguments ------------------------------------
731 !scalars
732  integer,intent(in) :: istwf_k,nx,ny,nz,ldx,ldy,ldz,ndat,npw
733 !arrays
734  integer,intent(in) :: kg_k(3,npw)
735  complex(spc),intent(in) :: iarrsph(npw*ndat)
736  complex(spc),intent(out) :: oarrbox(ldx*ldy*ldz*ndat)
737 
738 !Local variables-------------------------------
739 !scalars
740  integer,parameter :: me_g0=1
741  integer :: ix,ixinv,iy,iyinv,iz,izinv,dat,ipw,npwmin,pad_box,pad_sph,ifft,ifft_inv,ldxyz
742  !character(len=500) :: msg
743 !arrays
744  integer,allocatable :: ixinver(:),iyinver(:),izinver(:)
745 
746 ! *************************************************************************
747 
748 !In the case of special k-points, invariant under time-reversal,
749 !but not Gamma, initialize the inverse coordinates
750 !Remember indeed that
751 !u_k(G) = u_{k+G0}(G-G0); u_{-k}(G) = u_k(G)^*
752 !and therefore:
753 !u_{G0/2}(G) = u_{G0/2}(-G-G0)^*.
754  if (istwf_k>=2) then
755    ABI_MALLOC(ixinver,(nx))
756    ABI_MALLOC(iyinver,(ny))
757    ABI_MALLOC(izinver,(nz))
758    if ( ANY(istwf_k==(/2,4,6,8/)) ) then
759      ixinver(1)=1
760      do ix=2,nx
761        ixinver(ix)=nx+2-ix
762      end do
763    else
764      do ix=1,nx
765        ixinver(ix)=nx+1-ix
766      end do
767    end if
768    if (istwf_k>=2 .and. istwf_k<=5) then
769      iyinver(1)=1
770      do iy=2,ny
771        iyinver(iy)=ny+2-iy
772      end do
773    else
774      do iy=1,ny
775        iyinver(iy)=ny+1-iy
776      end do
777    end if
778    if ( ANY(istwf_k==(/2,3,6,7/)) ) then
779      izinver(1)=1
780      do iz=2,nz
781        izinver(iz)=nz+2-iz
782      end do
783    else
784      do iz=1,nz
785        izinver(iz)=nz+1-iz
786      end do
787    end if
788  end if
789 
790  ldxyz = ldx*ldy*ldz
791 
792  if (istwf_k==1) then
793 
794 !$OMP PARALLEL DO PRIVATE(pad_sph,pad_box,ix,iy,iz,ifft)
795    do dat=1,ndat
796      pad_sph = (dat-1)*npw
797      pad_box = (dat-1)*ldxyz
798      oarrbox(1+pad_box:ldxyz+pad_box) = czero_spc ! zero the sub-array
799      do ipw=1,npw
800        ix=kg_k(1,ipw); if (ix<0) ix=ix+nx; ix=ix+1
801        iy=kg_k(2,ipw); if (iy<0) iy=iy+ny; iy=iy+1
802        iz=kg_k(3,ipw); if (iz<0) iz=iz+nz; iz=iz+1
803        ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
804 #if defined __INTEL_COMPILER && defined HAVE_OPENMP
805        if (ifft==0) then
806          MSG_ERROR("prevent ifort+OMP from miscompiling this section on cronos")
807        end if
808 #endif
809        oarrbox(ifft+pad_box) = iarrsph(ipw+pad_sph)
810      end do
811    end do
812 
813  else if (istwf_k>=2) then
814    !
815    npwmin=1
816    if(istwf_k==2 .and. me_g0==1) then ! If gamma point, then oarrbox must be completed
817      do dat=1,ndat
818        pad_sph = (dat-1)*npw
819        pad_box = (dat-1)*ldxyz
820        oarrbox(1+pad_box) = REAL(iarrsph(1+pad_sph))
821      end do
822      npwmin=2
823    end if
824 
825 !$OMP PARALLEL DO PRIVATE(pad_sph,pad_box,ix,ixinv,iy,iyinv,iz,izinv,ifft)
826    do dat=1,ndat
827      pad_sph = (dat-1)*npw
828      pad_box = (dat-1)*ldxyz
829      oarrbox(npwmin+pad_box:ldxyz+pad_box) = czero_spc
830      do ipw=npwmin,npw
831        ix=kg_k(1,ipw); if(ix<0)ix=ix+nx; ix=ix+1
832        iy=kg_k(2,ipw); if(iy<0)iy=iy+ny; iy=iy+1
833        iz=kg_k(3,ipw); if(iz<0)iz=iz+nz; iz=iz+1
834        ifft = ix + (iy-1)*ldx + (iz-1)*ldx*ldy
835        ! Construct the coordinates of -k-G
836        ixinv=ixinver(ix); iyinv=iyinver(iy); izinv=izinver(iz)
837        ifft_inv = ixinv + (iyinv-1)*ldx + (izinv-1)*ldx*ldy
838 #if defined __INTEL_COMPILER && defined HAVE_OPENMP
839        if (ifft==0 .or. ifft_inv==0) then
840          MSG_ERROR("prevent ifort+OMP from miscompiling this section on cronos")
841        end if
842 #endif
843        oarrbox(ifft    +pad_box) =       iarrsph(ipw+pad_sph)
844        oarrbox(ifft_inv+pad_box) = CONJG(iarrsph(ipw+pad_sph))
845      end do
846    end do
847    !
848  else
849    MSG_ERROR("Wrong istwfk")
850  end if
851 
852  if (istwf_k>=2) then
853    ABI_FREE(ixinver)
854    ABI_FREE(iyinver)
855    ABI_FREE(izinver)
856  end if
857 
858 end subroutine cplx_gsph2box_spc

m_cplxtools/cplx_real_zdotc [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_real_zdotc

FUNCTION

   Perform a vector-vector operation defined as res = REAL (\Sigma (conjg(x)*y)) where x and y are n-element vectors.

INPUTS

  n = Specifies the number of elements in vector x and y
  x,y = Input arrays.

OUTPUT

  res=Real part of the scalar product.

PARENTS

CHILDREN

SOURCE

220 function cplx_real_zdotc(n,x,y) result(res)
221 
222 
223 !This section has been created automatically by the script Abilint (TD).
224 !Do not modify the following lines by hand.
225 #undef ABI_FUNC
226 #define ABI_FUNC 'cplx_real_zdotc'
227 !End of the abilint section
228 
229  implicit none
230 
231 !Arguments ------------------------------------
232 !scalars
233  integer,intent(in) :: n
234 !arrays
235  complex(dpc),intent(in) :: x(n)
236  complex(dpc),intent(in) :: y(n)
237  real(dp) :: res
238 
239 !Local variables-------------------------------
240  real(dp),external :: ddot
241 
242 ! *************************************************************************
243 
244  res = ddot(2*n,x,1,y,1)
245 
246 end function cplx_real_zdotc

m_cplxtools/cplx_setaug_zero_dpc [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_setaug_zero_dpc

FUNCTION

  Set to zero all elements of the array that are not in the FFT box.

INPUTS

 nx,ny,nz=physical dimensions of the FFT box
 ldx,ldy,ldx=memory dimension of arr
 ndat=number of FFTs

 SIDE EFFECT
  arr(ldx,ldy,ldz*ndat)= all entries in the augmented region are set to zero

PARENTS

      m_fft

CHILDREN

SOURCE

1134 subroutine cplx_setaug_zero_dpc(nx,ny,nz,ldx,ldy,ldz,ndat,arr)
1135 
1136 
1137 !This section has been created automatically by the script Abilint (TD).
1138 !Do not modify the following lines by hand.
1139 #undef ABI_FUNC
1140 #define ABI_FUNC 'cplx_setaug_zero_dpc'
1141 !End of the abilint section
1142 
1143  implicit none
1144 
1145 !Arguments ------------------------------------
1146 !scalars
1147  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat
1148 !arrays
1149  complex(dpc),intent(inout) :: arr(ldx,ldy,ldz*ndat)
1150 
1151 !Local variables-------------------------------
1152  integer :: iy,iz,dat,padat
1153 
1154 ! *************************************************************************
1155 
1156  if (nx /= ldx) then
1157    do iz=1,ldz*ndat
1158      do iy=1,ldy
1159        arr(nx+1:ldx,iy,iz) = czero_dpc
1160      end do
1161    end do
1162  end if
1163 
1164  if (ny /= ldy) then
1165    do iz=1,ldz*ndat
1166      arr(:,ny+1:ldy,iz) = czero_dpc
1167    end do
1168  end if
1169 
1170  if (nz /= ldz) then
1171    do dat=1,ndat
1172      padat = ldz*(dat-1)
1173      do iz=nz+1,ldz
1174        arr(:,:,iz+padat) = czero_dpc
1175      end do
1176    end do
1177  end if
1178 
1179 end subroutine cplx_setaug_zero_dpc

m_cplxtools/cplx_setaug_zero_spc [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_setaug_zero_spc

FUNCTION

  Set to zero all elements of the array that are not in the FFT box.

INPUTS

 nx,ny,nz=physical dimensions of the FFT box
 ldx,ldy,ldx=memory dimension of arr
 ndat=number of FFTs

 SIDE EFFECT
  arr(ldx,ldy,ldz*ndat)= all entries in the augmented region are set to zero

PARENTS

      m_fft

CHILDREN

SOURCE

1062 subroutine cplx_setaug_zero_spc(nx,ny,nz,ldx,ldy,ldz,ndat,arr)
1063 
1064 
1065 !This section has been created automatically by the script Abilint (TD).
1066 !Do not modify the following lines by hand.
1067 #undef ABI_FUNC
1068 #define ABI_FUNC 'cplx_setaug_zero_spc'
1069 !End of the abilint section
1070 
1071  implicit none
1072 
1073 !Arguments ------------------------------------
1074 !scalars
1075  integer,intent(in) :: nx,ny,nz,ldx,ldy,ldz,ndat
1076 !arrays
1077  complex(spc),intent(inout) :: arr(ldx,ldy,ldz*ndat)
1078 
1079 !Local variables-------------------------------
1080  integer :: iy,iz,dat,padat
1081 
1082 ! *************************************************************************
1083 
1084  if (nx /= ldx) then
1085    do iz=1,ldz*ndat
1086      do iy=1,ldy
1087        arr(nx+1:ldx,iy,iz) = czero_spc
1088      end do
1089    end do
1090  end if
1091 
1092  if (ny /= ldy) then
1093    do iz=1,ldz*ndat
1094      arr(:,ny+1:ldy,iz) = czero_spc
1095    end do
1096  end if
1097 
1098  if (nz /= ldz) then
1099    do dat=1,ndat
1100      padat = ldz*(dat-1)
1101      do iz=nz+1,ldz
1102        arr(:,:,iz+padat) = czero_spc
1103      end do
1104    end do
1105  end if
1106 
1107 end subroutine cplx_setaug_zero_spc

m_cplxtools/cplx_zaxpby [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_zaxpby

FUNCTION

  Scales two vectors, adds them to one another and stores result in the vector.
  y := a*x + b*y

INPUTS

 n = the number of elements in vectors x and y.
 a = Specifies the scalar a.
 x = Array.
 b = Specifies the scalar b.
 y = Array

OUTPUT

 y Contains the updated vector y.

PARENTS

CHILDREN

SOURCE

275 subroutine cplx_zaxpby(n,a,x,b,y)
276 
277 
278 !This section has been created automatically by the script Abilint (TD).
279 !Do not modify the following lines by hand.
280 #undef ABI_FUNC
281 #define ABI_FUNC 'cplx_zaxpby'
282 !End of the abilint section
283 
284  implicit none
285 
286 !Arguments ------------------------------------
287 !scalars
288  integer,intent(in) :: n
289  complex(dpc),intent(in) :: a,b
290 !arrays
291  complex(dpc),intent(in) :: x(n)
292  complex(dpc),intent(inout) :: y(n)
293 
294 ! *************************************************************************
295 
296 #ifdef HAVE_LINALG_AXPBY
297  call zaxpby(n, a, x, 1, b, y, 1)
298 #else
299  call zscal(n, b, y, 1)
300  call zaxpy(n, a, x, 1, y,1)
301 #endif
302 
303 end subroutine cplx_zaxpby

m_cplxtools/cplx_zgemm [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_zgemm

FUNCTION

  The ?gemm routines perform a matrix-matrix operation with general matrices.
  The operation is defined as C := alpha*op(A)*op(B) + beta*C,
  where:

  op(x) is one of op(x) = x, or op(x) = x', or op(x) = conjg(x'),

  alpha and beta are scalars,
  A, B and C are matrices:
  op(A) is an m-by-k matrix,
  op(B) is a k-by-n matrix,
  C is an m-by-n matrix.

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

413 subroutine cplx_zgemm(transa,transb,npws,ncola,ncolb,amat,bmat,cmat,alpha,beta)
414 
415 
416 !This section has been created automatically by the script Abilint (TD).
417 !Do not modify the following lines by hand.
418 #undef ABI_FUNC
419 #define ABI_FUNC 'cplx_zgemm'
420 !End of the abilint section
421 
422  implicit none
423 
424 !Arguments ------------------------------------
425 !scalars
426  integer,intent(in) :: npws,ncola,ncolb
427  complex(dpc),optional,intent(in) :: alpha,beta
428  character(len=1),intent(in) :: transa,transb
429 !arrays
430  complex(dpc),intent(in) :: amat(npws*ncola)
431  complex(dpc),intent(in) :: bmat(npws*ncolb)
432  complex(dpc),intent(inout) :: cmat(*)
433 
434 !Local variables-------------------------------
435 !scalars
436  integer :: mm,nn,kk,lda,ldb,ldc
437  complex(dpc) :: my_alpha,my_beta
438 
439 ! *************************************************************************
440 
441  lda = npws
442  ldb = npws
443 
444  mm  = npws
445  nn  = ncolb
446  kk  = ncola
447 
448  if (toupper(transa) /= 'N') then
449    mm = ncola
450    kk = npws
451  end if
452  if (toupper(transb) /= 'N') nn = npws
453 
454  ldc = mm
455 
456  my_alpha = cone_dpc;  if (PRESENT(alpha)) my_alpha = alpha
457  my_beta  = czero_dpc; if (PRESENT(beta))  my_beta  = beta
458 
459  call ZGEMM(transa,transb,mm,nn,kk,my_alpha,amat,lda,bmat,ldb,my_beta,cmat,ldc)
460 
461 end subroutine cplx_zgemm

m_cplxtools/cplx_zgemv [ Functions ]

[ Top ] [ m_cplxtools ] [ Functions ]

NAME

  cplx_zgemv

FUNCTION

 The ?gemv routines perform a matrix-vector operation defined as

 y := alpha*A*x + beta*y,
 or
 y := alpha*A'*x + beta*y,
 or
 y := alpha*conjg(A')*x + beta*y,

 where: alpha and beta are scalars, x and y are vectors, A is an m-by-n matrix.

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

333 subroutine cplx_zgemv(trans,nrows,ncols,mat,vec,matvec,alpha,beta)
334 
335 
336 !This section has been created automatically by the script Abilint (TD).
337 !Do not modify the following lines by hand.
338 #undef ABI_FUNC
339 #define ABI_FUNC 'cplx_zgemv'
340 !End of the abilint section
341 
342  implicit none
343 
344 !Arguments ------------------------------------
345 !scalars
346  integer,intent(in) :: nrows,ncols
347  complex(dpc),optional,intent(in) :: alpha,beta
348  character(len=1),intent(in) :: trans
349 !arrays
350  complex(dpc),intent(in) :: mat(nrows*ncols)
351  complex(dpc),intent(in) :: vec(*)
352  complex(dpc),intent(inout) :: matvec(*)
353 
354 !Local variables-------------------------------
355 !scalars
356  integer :: mm,nn,kk,lda,ldb,ldc
357  complex(dpc) :: my_alpha,my_beta
358 
359 ! *************************************************************************
360 
361  lda = nrows
362  mm  = nrows
363  nn  = 1
364  kk  = ncols
365 
366  if (toupper(trans) /= 'N') then
367    mm = ncols
368    kk = nrows
369  end if
370 
371  ldb = kk
372  ldc = mm
373 
374  my_alpha = cone_dpc;  if (PRESENT(alpha)) my_alpha = alpha
375  my_beta  = czero_dpc; if (PRESENT(beta))  my_beta  = beta
376 
377  call ZGEMM(trans,"N",mm,nn,kk,my_alpha,mat,lda,vec,ldb,my_beta,matvec,ldc)
378 
379  !call ZGEMV(trans,mm,nn,my_alpha,mat,lda,vec,1,my_beta,matvec,1)
380 
381 end subroutine cplx_zgemv