1!
2! Copyright (C) 2002-2011 Quantum ESPRESSO groups
3! This file is distributed under the terms of the
4! GNU General Public License. See the file `License'
5! in the root directory of the present distribution,
6! or http://www.gnu.org/copyleft/gpl.txt .
7!
8!-----------------------------------------------------------------------
9      SUBROUTINE initbox ( tau0, alat, at, ainv, taub, irb )
10!-----------------------------------------------------------------------
11!
12!     sets the indexes irb and positions taub for the small boxes
13!     around atoms
14!
15      USE kinds,                    ONLY: DP
16      USE ions_base,                ONLY: nat
17      USE control_flags,            ONLY: iverbosity
18      USE io_global,                ONLY: stdout
19      USE mp_global,                ONLY: nproc_bgrp, me_bgrp, intra_bgrp_comm
20      USE fft_base,                 ONLY: dfftb, dfftp, fft_type_descriptor
21      USE fft_smallbox_type,        ONLY: fft_box_set
22
23      IMPLICIT NONE
24! input
25      REAL(DP), INTENT(in)  :: tau0(3,nat), at(3,3), ainv(3,3), alat
26! output
27      INTEGER,  INTENT(out) :: irb(3,nat)
28      REAL(DP), INTENT(out) :: taub(3,nat)
29! local
30      REAL(DP) :: x(3), xmod
31      INTEGER  :: nr(3), nrb(3), xint, ia, i
32!
33      IF ( dfftb%nr1 < 1) CALL errore ('initbox', 'incorrect value for box grid dimensions', 1)
34      IF ( dfftb%nr2 < 1) CALL errore ('initbox', 'incorrect value for box grid dimensions', 2)
35      IF ( dfftb%nr3 < 1) CALL errore ('initbox', 'incorrect value for box grid dimensions', 3)
36
37      nr (1)=dfftp%nr1
38      nr (2)=dfftp%nr2
39      nr (3)=dfftp%nr3
40      nrb(1)=dfftb%nr1
41      nrb(2)=dfftb%nr2
42      nrb(3)=dfftb%nr3
43!
44      DO ia=1,nat
45!
46            DO i=1,3
47!
48! bring atomic positions to crystal axis
49!
50               x(i) = ainv(i,1)*tau0(1,ia) +                         &
51     &                ainv(i,2)*tau0(2,ia) +                         &
52     &                ainv(i,3)*tau0(3,ia)
53!
54! bring x in the range between 0 and 1
55!
56               x(i) = MOD(x(i),1.d0)
57               IF (x(i).LT.0.d0) x(i)=x(i)+1.d0
58!
59! case of nrb(i) even
60!
61               IF (MOD(nrb(i),2).EQ.0) THEN
62!
63! find irb = index of the grid point at the corner of the small box
64!           (the indices of the small box run from irb to irb+nrb-1)
65!
66                  xint=INT(x(i)*nr(i))
67                  irb (i,ia)=xint+1-nrb(i)/2+1
68                  IF(irb(i,ia).LT.1) irb(i,ia)=irb(i,ia)+nr(i)
69!
70! x(i) are the atomic positions in crystal coordinates, where the
71! "crystal lattice" is the small box lattice and the origin is at
72! the corner of the small box. Used to calculate phases exp(iG*taub)
73!
74                  xmod=x(i)*nr(i)-xint
75                  x(i)=(xmod+nrb(i)/2-1)/nr(i)
76               ELSE
77!
78! case of nrb(i) odd - see above for comments
79!
80                  xint=NINT(x(i)*nr(i))
81                  irb (i,ia)=xint+1-(nrb(i)-1)/2
82                  IF(irb(i,ia).LT.1) irb(i,ia)=irb(i,ia)+nr(i)
83                  xmod=x(i)*nr(i)-xint
84                  x(i)=(xmod+(nrb(i)-1)/2)/nr(i)
85               END IF
86            END DO
87!
88! bring back taub in cartesian coordinates
89!
90            DO i=1,3
91               taub(i,ia)=(x(1)*at(i,1) + x(2)*at(i,2) + x(3)*at(i,3))*alat
92            END DO
93      END DO
94
95      ! initialize FFT descriptor
96
97      CALL fft_box_set( dfftb, nat, irb, dfftp )
98
99      IF( iverbosity > 1 ) THEN
100           DO ia=1,nat
101              WRITE( stdout,2000) ia, (irb(i,ia),i=1,3)
102           END DO
1032000       FORMAT(2x, 'atom= ', i3, ' irb1= ', i3, ' irb2= ', i3, ' irb3= ', i3)
104      ENDIF
105
106#if defined(__MPI)
107      !
108      ! for processor that do not call fft on the box
109      ! artificially start the clock
110      !
111      CALL start_clock( 'fftb' )
112      CALL stop_clock( 'fftb' )
113      !
114#endif
115!
116      RETURN
117   END SUBROUTINE initbox
118!
119!-----------------------------------------------------------------------
120      SUBROUTINE phbox( taub, iverbosity, eigrb )
121!-----------------------------------------------------------------------
122!     calculates the phase factors for the g's of the little box
123!     eigrt=exp(-i*g*tau) .
124!     Uses the same logic for fast calculation as in phfac
125!
126      USE kinds,         only: DP
127      use io_global,     only: stdout
128      use ions_base,     only: nsp, na, nat
129      use cell_base,     only: r_to_s
130      use cp_interfaces, only: phfacs
131      use small_box,     only: bgb, alatb
132      use smallbox_gvec, only: ngb, mill_b
133      use fft_base,      only: dfftb
134!
135      IMPLICIT NONE
136      REAL(DP), INTENT(IN)    :: taub(3,nat)
137      COMPLEX(DP), INTENT(OUT) :: eigrb(ngb,nat)
138      INTEGER, INTENT(IN) :: iverbosity
139! local
140      REAL(DP)    :: ainvb(3,3)
141      integer :: i, ia, ig
142      complex(dp), allocatable:: ei1b(:,:), ei2b(:,:), ei3b(:,:)
143      real(dp), allocatable :: taus(:,:)
144!
145      allocate(ei1b(-dfftb%nr1:dfftb%nr1,nat))
146      allocate(ei2b(-dfftb%nr2:dfftb%nr2,nat))
147      allocate(ei3b(-dfftb%nr3:dfftb%nr3,nat))
148      allocate( taus( 3, nat ) )
149!
150      if(iverbosity > 2) then
151         WRITE( stdout,*) ' phbox: taub '
152         WRITE( stdout,*) ( (taub(i,ia), i=1, 3 ), ia=1, nat )
153      endif
154
155      ainvb(1,:) = bgb(:,1)/alatb
156      ainvb(2,:) = bgb(:,2)/alatb
157      ainvb(3,:) = bgb(:,3)/alatb
158
159      CALL r_to_s( taub, taus, nat, ainvb )
160      CALL phfacs( ei1b, ei2b, ei3b, eigrb, mill_b, taus, dfftb%nr1,dfftb%nr2,dfftb%nr3, nat )
161!
162      if(iverbosity > 2) then
163         WRITE( stdout,*)
164         do ia=1,nat
165            WRITE( stdout,'(33x,a,i4)') ' ei1b, ei2b, ei3b (ia)',ia
166            do ig=1,4
167               WRITE( stdout,'(6f9.4)') ei1b(ig,ia),ei2b(ig,ia),ei3b(ig,ia)
168            end do
169            WRITE( stdout,*)
170         end do
171      endif
172!
173      deallocate(ei3b)
174      deallocate(ei2b)
175      deallocate(ei1b)
176      deallocate( taus )
177!
178      RETURN
179      END SUBROUTINE phbox
180
181