TABLE OF CONTENTS


ABINIT/phase [ Functions ]

[ Top ] [ Functions ]

NAME

 phase

FUNCTION

 Compute ph(ig)=$\exp(\pi\ i \ n/ngfft)$ for n=0,...,ngfft/2,-ngfft/2+1,...,-1
 while ig runs from 1 to ngfft.

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR)
 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 .

INPUTS

  ngfft=number of points

OUTPUT

  ph(2*ngfft)=phase array (complex)

NOTES

 XG 990504 : changed the formulation, in order to preserve
 the invariance between n and -n, that was broken for n=ngfft/2 if ngfft even.
 Simply suppresses the corresponding sine.

PARENTS

      xcden,xcpot

CHILDREN

SOURCE

34 #if defined HAVE_CONFIG_H
35 #include "config.h"
36 #endif
37 
38 #include "abi_common.h"
39 
40 
41 subroutine phase(ngfft,ph)
42 
43  use defs_basis
44  use m_profiling_abi
45 
46 !This section has been created automatically by the script Abilint (TD).
47 !Do not modify the following lines by hand.
48 #undef ABI_FUNC
49 #define ABI_FUNC 'phase'
50 !End of the abilint section
51 
52  implicit none
53 
54 !Arguments ------------------------------------
55 !scalars
56  integer,intent(in) :: ngfft
57 !arrays
58  real(dp),intent(out) :: ph(2*ngfft)
59 
60 !Local variables-------------------------------
61 !scalars
62  integer :: id,ig,nn
63  real(dp) :: arg,fac
64 
65 ! *************************************************************************
66 
67  id=ngfft/2+2
68  fac=pi/dble(ngfft)
69  do ig=1,ngfft
70    nn=ig-1-(ig/id)*ngfft
71    arg=fac*dble(nn)
72    ph(2*ig-1)=cos(arg)
73    ph(2*ig)  =sin(arg)
74 
75  end do
76 !XG 990504 Here zero the corresponding sine
77  if((ngfft/2)*2==ngfft) ph(2*(id-1))=zero
78 
79 end subroutine phase