TABLE OF CONTENTS


ABINIT/scalewf_nonlop [ Functions ]

[ Top ] [ Functions ]

NAME

 scalewf_nonlop

FUNCTION

 At the start of nonlop (or similar routines), as well as its end,
 the wavefunctions, when stored with istwfk/=2,
 need to be scaled (by a factor of 2 or 1/2),
 except for the G=0 component.
 Only the first spinor component is to be modified.

COPYRIGHT

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

INPUTS

  istwf_k=storage mode of the vector
  mpi_enreg=informations about MPI parallelization
  npw=number of planewaves
  option=1 multiply by 2
        =2 multiply by 1/2

OUTPUT

  (see side effects)

SIDE EFFECTS

  vect(2,npw)=vector that is rescaled

NOTES

  XG030513 : MPIWF One should pay attention to the
  G=0 component, that will be only one one proc...

PARENTS

      nonlop_pl

CHILDREN

SOURCE

 43 #if defined HAVE_CONFIG_H
 44 #include "config.h"
 45 #endif
 46 
 47 #include "abi_common.h"
 48 
 49 
 50 subroutine scalewf_nonlop(istwf_k,mpi_enreg,npw,option,vect)
 51 
 52  use defs_basis
 53  use defs_abitypes
 54  use m_errors
 55  use m_profiling_abi
 56 
 57 !This section has been created automatically by the script Abilint (TD).
 58 !Do not modify the following lines by hand.
 59 #undef ABI_FUNC
 60 #define ABI_FUNC 'scalewf_nonlop'
 61 !End of the abilint section
 62 
 63  implicit none
 64 
 65 !Arguments ------------------------------------
 66 !This type is defined in defs_mpi
 67 !scalars
 68  integer,intent(in) :: istwf_k,npw,option
 69  type(MPI_type),intent(in) :: mpi_enreg
 70 !arrays
 71  real(dp),intent(inout) :: vect(2,npw)
 72 
 73 !Local variables-------------------------------
 74 !scalars
 75  integer :: ipw
 76  real(dp) :: scale
 77  character(len=500) :: message
 78 
 79 ! *************************************************************************
 80 
 81  DBG_ENTER("COLL")
 82 
 83  if(istwf_k/=1)then
 84 
 85    if(option/=1 .and. option/=2)then
 86      write(message,'(a,a,a,i0)')&
 87 &     'The argument option should be 1 or 2,',ch10,&
 88 &     'however, option=',option
 89      MSG_BUG(message)
 90    end if
 91 
 92    scale=two
 93    if(option==2)scale=half
 94 
 95 !  Storage for the Gamma point. The component of the G=0 vector
 96 !  should not be scaled, and no G=0 imaginary part is allowed.
 97    if(istwf_k==2)then
 98      if (mpi_enreg%me_g0==1) then
 99        vect(2,1)=zero
100 !$OMP PARALLEL DO 
101        do ipw=2,npw
102          vect(1,ipw)=scale*vect(1,ipw)
103          vect(2,ipw)=scale*vect(2,ipw)
104        end do
105 !$OMP END PARALLEL DO 
106      else
107 !$OMP PARALLEL DO 
108        do ipw=1,npw
109          vect(1,ipw)=scale*vect(1,ipw)
110          vect(2,ipw)=scale*vect(2,ipw)
111        end do
112 !$OMP END PARALLEL DO 
113      end if
114    end if
115 
116 !  Other storage modes, for k points with time-reversal symmetry.
117 !  All components should be scaled.
118    if(istwf_k>2)then
119 !$OMP PARALLEL DO 
120      do ipw=1,npw
121        vect(1,ipw)=scale*vect(1,ipw)
122        vect(2,ipw)=scale*vect(2,ipw)
123      end do
124 !$OMP END PARALLEL DO 
125    end if
126 
127  end if ! istwf_k/=1
128 
129  DBG_EXIT("COLL")
130 
131 end subroutine scalewf_nonlop