1!
2! Copyright (C) 2001-2020 Quantum ESPRESSO Foundation
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!-----------------------------------------------------------------------
9subroutine star_q1(xq, at, bg, nsym, s, invs, nq, sxq, isq, imq, verbosity,&
10     t_rev )
11  !-----------------------------------------------------------------------
12  ! FIXME: MERGE WITH STAR_Q
13  ! generate the star of q vectors that are equivalent to the input one
14  ! NB: input s(:,:,1:nsym) must contain all crystal symmetries,
15  ! i.e. not those of the small-qroup of q only
16  !
17  USE io_global,  ONLY : stdout
18  USE kinds, only : DP
19  implicit none
20  !
21  real(DP), parameter :: accep=1.e-5_dp
22
23  integer, intent(in) :: nsym, s (3, 3, 48), invs(48)
24  integer, intent(in) :: t_rev(48)
25  ! nsym matrices of symmetry operations
26  ! invs: list of inverse operation indices
27  ! t_rev: if a given simmetry operation needs time reversal
28  real(DP), intent(in) :: xq (3), at (3, 3), bg (3, 3)
29  ! xq: q vector
30  ! at: direct lattice vectors
31  ! bg: reciprocal lattice vectors
32  !
33  integer, intent(out) :: nq, isq (48), imq
34  ! nq  : degeneracy of the star of q
35  ! isq : index of q in the star for a given sym
36  ! imq : index of -q in the star (0 if not present)
37
38  real(DP), intent(out) :: sxq (3, 48)
39  ! list of vectors in the star of q
40  logical, intent(in) :: verbosity
41  ! if true prints several messages.
42  !
43  integer :: nsq (48), isym, ism1, iq, i
44  ! number of symmetry ops. of bravais lattice
45  ! counters on symmetry ops.
46  ! index of inverse of isym
47  ! counters
48  real(DP) :: saq (3, 48), aq (3), raq (3), zero (3)
49  ! auxiliary list of q (crystal coordinates)
50  ! input q in crystal coordinates
51  ! rotated q in crystal coordinates
52  ! coordinates of fractionary translations
53  ! a zero vector: used in eqvect
54
55  logical, external :: eqvect
56  ! function used to compare two vectors
57  !
58  !
59  zero(:) = 0.d0
60  saq(:, :) = 0.0d0
61  !
62  ! go to  crystal coordinates
63  !
64  do i = 1, 3
65     aq(i) = xq(1) * at(1,i) + xq(2) * at(2,i) + xq(3) * at(3,i)
66  enddo
67  !
68  ! create the list of rotated q
69  !
70  do i = 1, 48
71     nsq (i) = 0
72     isq (i) = 0
73  enddo
74  nq = 0
75  do isym = 1, nsym
76     ism1 = invs (isym)
77     do i = 1, 3
78        raq (i) = s (i, 1, ism1) * aq (1) &
79                + s (i, 2, ism1) * aq (2) &
80                + s (i, 3, ism1) * aq (3)
81     enddo
82     IF (t_rev(isym)==1) raq = -raq
83     do i = 1, 3
84        sxq (i, 48) = bg (i, 1) * raq (1) &
85                    + bg (i, 2) * raq (2) &
86                    + bg (i, 3) * raq (3)
87     enddo
88     do iq = 1, nq
89        if (eqvect (raq, saq (1, iq), zero, accep) ) then
90           isq (isym) = iq
91           nsq (iq) = nsq (iq) + 1
92        endif
93     enddo
94     if (isq (isym) == 0) then
95        nq = nq + 1
96        nsq (nq) = 1
97        isq (isym) = nq
98        saq(:,nq) = raq(:)
99        do i = 1, 3
100           sxq (i, nq) = bg (i, 1) * saq (1, nq) &
101                       + bg (i, 2) * saq (2, nq) &
102                       + bg (i, 3) * saq (3, nq)
103        enddo
104     endif
105  enddo
106  !
107  ! set imq index if needed and check star degeneracy
108  !
109  raq (:) = - aq(:)
110  imq = 0
111  do iq = 1, nq
112     if (eqvect (raq, saq (1, iq), zero, accep) ) imq = iq
113     if (nsq(iq)*nq /= nsym) call errore ('star_q', 'wrong degeneracy', iq)
114  enddo
115  !
116  ! writes star of q
117  !
118  IF (verbosity) THEN
119  WRITE( stdout, * )
120  WRITE( stdout, '(5x,a,i4)') 'Number of q in the star = ', nq
121  WRITE( stdout, '(5x,a)') 'List of q in the star:'
122  WRITE( stdout, '(7x,i4,3f14.9)') (iq, (sxq(i,iq), i=1,3), iq=1,nq)
123  if (imq == 0) then
124     WRITE( stdout, '(5x,a)') 'In addition there is the -q list: '
125     WRITE( stdout, '(7x,i4,3f14.9)') (iq, (-sxq(i,iq), i=1,3), iq=1,nq)
126  endif
127  ENDIF
128  return
129end subroutine star_q1
130!-----------------------------------------------------------------------
131subroutine star_q (xq, at, bg, nsym, s, invs, nq, sxq, isq, imq, verbosity)
132  !-----------------------------------------------------------------------
133  ! generate the star of q vectors that are equivalent to the input one
134  ! NB: input s(:,:,1:nsym) must contain all crystal symmetries,
135  ! i.e. not those of the small-qroup of q only
136  !
137  USE io_global,  ONLY : stdout
138  USE kinds, only : DP
139  implicit none
140  !
141  real(DP), parameter :: accep=1.e-5_dp
142
143  integer, intent(in) :: nsym, s (3, 3, 48), invs(48)
144  ! nsym matrices of symmetry operations
145  ! invs: list of inverse operation indices
146  real(DP), intent(in) :: xq (3), at (3, 3), bg (3, 3)
147  ! xq: q vector
148  ! at: direct lattice vectors
149  ! bg: reciprocal lattice vectors
150  !
151  integer, intent(out) :: nq, isq (48), imq
152  ! nq  : degeneracy of the star of q
153  ! isq : index of q in the star for a given sym
154  ! imq : index of -q in the star (0 if not present)
155
156  real(DP), intent(out) :: sxq (3, 48)
157  ! list of vectors in the star of q
158  logical, intent(in) :: verbosity
159  ! if true prints several messages.
160  !
161  integer :: nsq (48), isym, ism1, iq, i, t_rev(48)
162  ! number of symmetry ops. of bravais lattice
163  ! counters on symmetry ops.
164  ! index of inverse of isym
165  ! counters
166  ! t_rev variable, says if a given simmetry operation
167  ! needs the time reversal to be a symmetry of the crystal.
168  real(DP) :: saq (3, 48), aq (3), raq (3), zero (3)
169  ! auxiliary list of q (crystal coordinates)
170  ! input q in crystal coordinates
171  ! rotated q in crystal coordinates
172  ! coordinates of fractionary translations
173  ! a zero vector: used in eqvect
174
175  logical, external :: eqvect
176  ! function used to compare two vectors
177  !
178  !
179  zero(:) = 0.d0
180  saq(:, :) = 0.0d0
181  !
182  ! go to  crystal coordinates
183  !
184  do i = 1, 3
185     aq(i) = xq(1) * at(1,i) + xq(2) * at(2,i) + xq(3) * at(3,i)
186  enddo
187  !
188  ! create the list of rotated q
189  !
190  do i = 1, 48
191     nsq (i) = 0
192     isq (i) = 0
193  enddo
194  nq = 0
195  do isym = 1, nsym
196     ism1 = invs (isym)
197     do i = 1, 3
198        raq (i) = s (i, 1, ism1) * aq (1) &
199                + s (i, 2, ism1) * aq (2) &
200                + s (i, 3, ism1) * aq (3)
201     enddo
202     do i = 1, 3
203        sxq (i, 48) = bg (i, 1) * raq (1) &
204                    + bg (i, 2) * raq (2) &
205                    + bg (i, 3) * raq (3)
206     enddo
207     do iq = 1, nq
208        if (eqvect (raq, saq (1, iq), zero, accep) ) then
209           isq (isym) = iq
210           nsq (iq) = nsq (iq) + 1
211        endif
212     enddo
213     if (isq (isym) == 0) then
214        nq = nq + 1
215        nsq (nq) = 1
216        isq (isym) = nq
217        saq(:,nq) = raq(:)
218        do i = 1, 3
219           sxq (i, nq) = bg (i, 1) * saq (1, nq) &
220                       + bg (i, 2) * saq (2, nq) &
221                       + bg (i, 3) * saq (3, nq)
222        enddo
223     endif
224  enddo
225  !
226  ! set imq index if needed and check star degeneracy
227  !
228  raq (:) = - aq(:)
229  imq = 0
230  do iq = 1, nq
231     if (eqvect (raq, saq (1, iq), zero, accep) ) imq = iq
232     if (nsq(iq)*nq /= nsym) call errore ('star_q', 'wrong degeneracy', iq)
233  enddo
234  !
235  ! writes star of q
236  !
237  IF (verbosity) THEN
238  WRITE( stdout, * )
239  WRITE( stdout, '(5x,a,i4)') 'Number of q in the star = ', nq
240  WRITE( stdout, '(5x,a)') 'List of q in the star:'
241  WRITE( stdout, '(7x,i4,3f14.9)') (iq, (sxq(i,iq), i=1,3), iq=1,nq)
242  if (imq == 0) then
243     WRITE( stdout, '(5x,a)') 'In addition there is the -q list: '
244     WRITE( stdout, '(7x,i4,3f14.9)') (iq, (-sxq(i,iq), i=1,3), iq=1,nq)
245  endif
246  ENDIF
247  return
248end subroutine star_q
249