Line data Source code
1 : !!****m* ABINIT/m_paw_sym
2 : !! NAME
3 : !! m_paw_sym
4 : !!
5 : !! FUNCTION
6 : !! This module contains several routines related to the use of symmetries in the PAW approach.
7 : !!
8 : !! COPYRIGHT
9 : !! Copyright (C) 2018-2026 ABINIT group (MG)
10 : !! This file is distributed under the terms of the
11 : !! GNU General Public License, see ~abinit/COPYING
12 : !! or http://www.gnu.org/copyleft/gpl.txt .
13 : !!
14 : !! SOURCE
15 :
16 : #if defined HAVE_CONFIG_H
17 : #include "config.h"
18 : #endif
19 :
20 : #include "abi_common.h"
21 :
22 : MODULE m_paw_sym
23 :
24 : use defs_basis
25 : use m_abicore
26 : use m_errors
27 :
28 : use m_crystal, only : crystal_t
29 : use m_pawang, only : pawang_type
30 : use m_pawtab, only : pawtab_type
31 : use m_pawcprj, only : pawcprj_type, pawcprj_alloc, pawcprj_free, pawcprj_copy
32 : use m_bz_mesh, only : kmesh_t
33 :
34 : implicit none
35 :
36 : private
37 :
38 : !public procedures.
39 : public :: paw_symcprj ! Symetrize the projections cprj=<n,k|p_i> (p_i=NL PAW projector) - in-place version
40 : public :: paw_symcprj_op ! Symetrize the projections cprj=<n,k|p_i> (p_i=NL PAW projector) - out-of-place version
41 :
42 : CONTAINS !========================================================================================
43 : !!***
44 :
45 : !----------------------------------------------------------------------
46 :
47 : !!****f* m_paw_sym/paw_symcprj
48 : !! NAME
49 : !! paw_symcprj
50 : !!
51 : !! FUNCTION
52 : !! Symetrize the projections cprj=<n,k|p_i> (p_i=NL PAW projector) - in-place version
53 : !!
54 : !! INPUTS
55 : !! ik_ibz=The index of the k-point in the full BZ where the matrix elements have to be symmetrized.
56 : !! nspinor=Number of spinorial components
57 : !! nband_k=Number of bands stored in cprjnk_kibz for this k-point.
58 : !! Cryst<crystal_t>=data type gathering information on unit cell and symmetries.
59 : !! %ntypat=number of type of atoms
60 : !! %natom=number of atoms in the unit cell
61 : !! %typat(natom)=type of each atom
62 : !! %indsym(4,nsym,natom)=indirect indexing array:
63 : !! for each isym,iatom, fourth element is label of atom into which iatom is sent by the INVERSE of the
64 : !! symmetry operation symrel(isym); first three elements are the primitive translations that must be subtracted
65 : !! after the transformation to get back to the original unit cell.
66 : !! Kmesh<kmesh_t>: datatype gathering information on the k-point sampling.
67 : !! %nbz=number of k-points in the full Brillouin zone
68 : !! %nibz=number of k-points in the irreducible wedge
69 : !! %tab(nkbz)=table giving for each k-point in the BZ (array kbz), the corresponding irred. point in the IBZ.
70 : !! i.e k_BZ = (IS) kIBZ where S is one of the symrec operations and I is the inversion or the identity
71 : !! %tabi(nkbz)=for each k-point in the BZ defines whether inversion has to be considered in the
72 : !! relation k_BZ=(IS) k_IBZ (1 => only S; -1 => -S)
73 : !! %tabo(nkbz)= the symmetry operation S that takes k_IBZ to each k_BZ
74 : !! Pawtab(Cryst%ntypat) <type(pawtab_type)>=paw tabulated starting data.
75 : !! Pawang <type(pawang_type)>=paw angular mesh and related data
76 : !! %lmax=Max angular momentum included in the PAW datasets used. mentioned at the second line of the psp file
77 : !! %zarot(2*lmax+1,2*lmax+1,lmax+1,nsym)=coefficients of the transformation of real spherical
78 : !! harmonics under the symmetry operations.
79 : !!
80 : !! NOTES
81 : !! Derivatives are not symmetrized.
82 : !!
83 : !! OUTPUT
84 : !!
85 : !! SOURCE
86 :
87 127533 : subroutine paw_symcprj(ik_bz,nspinor,nband_k,Cryst,Kmesh,Pawtab,Pawang,Cprj_bz)
88 :
89 : !Arguments ------------------------------------
90 : !scalars
91 : integer,intent(in) :: nspinor,nband_k,ik_bz
92 : type(crystal_t),intent(in) :: Cryst
93 : type(kmesh_t),intent(in) :: Kmesh
94 : type(Pawang_type),intent(in) :: Pawang
95 : !arrays
96 : type(Pawtab_type),target,intent(in) :: Pawtab(Cryst%ntypat)
97 : type(pawcprj_type),intent(inout) :: Cprj_bz(Cryst%natom,nspinor*nband_k)
98 :
99 : !Local variables-------------------------------
100 : !scalars
101 : integer :: iatom,iat_sym,iband,ibsp_bz
102 : integer :: ibsp_ibz,ik_ibz,indexj,ispinor,isym,itim
103 : integer :: itypat,jl,jl0,jlmn,jln,jln0,jlpm,jm,jn,lmax,mm,ncpgr
104 : real(dp) :: arg,wtk
105 : logical :: isirred
106 : !arrays
107 255066 : integer :: r0(3),nlmn_atom(Cryst%natom)
108 127533 : integer,ABI_CONTIGUOUS pointer :: indlmn(:,:)
109 23045 : real(dp) :: dum(2,nspinor),kbz(3),kirr(3),phase(2),swp(2),tmp(2,nspinor)
110 127533 : real(dp),allocatable :: DS_mmpl(:,:,:)
111 1828302 : type(pawcprj_type) :: Cprjnk_kibz(Cryst%natom,nspinor*nband_k)
112 :
113 : ! *********************************************************************
114 :
115 127533 : ncpgr = Cprj_bz(1,1)%ncpgr
116 127533 : ABI_CHECK(ncpgr==0,"Derivatives of cprj are not coded")
117 :
118 : !Get the index of the IBZ image associated to the BZ k-point ik_bz and related simmetry.
119 127533 : call Kmesh%get_BZ_item(ik_bz,kbz,ik_ibz,isym,itim,isirred=isirred)
120 :
121 281931 : if (isirred) RETURN ! It is a point in the IBZ, Symmetrization is not needed.
122 : !
123 : !The corresponding point kirr in the IBZ.
124 104488 : call kmesh%get_IBZ_item(ik_ibz,kirr,wtk)
125 :
126 : !Local copy.
127 316656 : do iatom=1,Cryst%natom
128 316656 : nlmn_atom(iatom)=Pawtab(Cryst%typat(iatom))%lmn_size
129 : end do
130 :
131 104488 : call pawcprj_alloc(Cprjnk_kibz,ncpgr,nlmn_atom)
132 104488 : call pawcprj_copy(Cprj_bz,Cprjnk_kibz)
133 : !
134 : !=== DS_mmpl is the rotation matrix for real spherical harmonics associated to symrec(:,:,isym) ===
135 : !* Note the convention used by Blanco in Eq. 27 : DS_mmp multiply spherical harmonics as row vectors
136 104488 : lmax=Pawang%l_max-1 ! l_max is Max l+1
137 626928 : ABI_MALLOC(DS_mmpl,(2*lmax+1,2*lmax+1,lmax+1))
138 2956216 : DS_mmpl=Pawang%zarot(:,:,:,isym)
139 : !
140 : !===========================================
141 : !==== Loop over atoms to be symmetrized ====
142 : !===========================================
143 316656 : do iatom=1,Cryst%natom
144 212168 : itypat=Cryst%typat(iatom)
145 212168 : iat_sym=Cryst%indsym(4,isym,iatom)
146 212168 : indlmn => Pawtab(itypat)%indlmn
147 848672 : r0=Cryst%indsym(1:3,isym,iatom) ! R^{-1} (xred(:,iatom)-tnons) = xred(:,iat_sym) + r0.
148 848672 : arg=two_pi*dot_product(kirr,r0)
149 212168 : phase(1)=COS(arg)
150 212168 : phase(2)=SIN(arg)
151 : !
152 : ! Loop over the (jl,jm,jn) components to be symmetrized.
153 212168 : jl0=-1; jln0=-1; indexj=1
154 2018950 : do jlmn=1,Pawtab(itypat)%lmn_size
155 1702294 : jl =indlmn(1,jlmn)
156 1702294 : jm =indlmn(2,jlmn)
157 1702294 : jn =indlmn(3,jlmn)
158 1702294 : jln =indlmn(5,jlmn)
159 1702294 : jlpm=1+jl+jm
160 1702294 : if (jln/=jln0) indexj=indexj+2*jl0+1
161 : !
162 : ! === For each band, calculate contribution due to rotated real spherical harmonics ===
163 : ! FIXME check this expression; according to Blanco I should have D(S^-1} but it seems D(S) is correct
164 : ! Recheck spinorial case, presently is wrong
165 1702294 : ibsp_ibz=0
166 1702294 : ibsp_bz=0
167 3404588 : do iband=1,nband_k
168 :
169 6809176 : tmp(:,:)=zero
170 3404588 : do ispinor=1,nspinor
171 1702294 : ibsp_ibz=ibsp_ibz+1
172 7676478 : do mm=1,2*jl+1
173 4271890 : tmp(1,ispinor)=tmp(1,ispinor)+DS_mmpl(mm,jlpm,jl+1)*Cprjnk_kibz(iat_sym,ibsp_ibz)%cp(1,indexj+mm)
174 5974184 : tmp(2,ispinor)=tmp(2,ispinor)+DS_mmpl(mm,jlpm,jl+1)*Cprjnk_kibz(iat_sym,ibsp_ibz)%cp(2,indexj+mm)
175 : end do
176 : end do !ispinor
177 : !
178 : ! * Apply the phase to account if the symmetric atom belongs to a different unit cell.
179 3404588 : do ispinor=1,nspinor
180 1702294 : dum(1,ispinor)=tmp(1,ispinor)*phase(1)-tmp(2,ispinor)*phase(2)
181 3404588 : dum(2,ispinor)=tmp(1,ispinor)*phase(2)+tmp(2,ispinor)*phase(1)
182 : end do
183 : !
184 : ! * If required, apply time-reversal symmetry to retrieve the correct point in the BZ.
185 1702294 : if (itim==2) then
186 0 : if (nspinor==1) then
187 0 : dum(2,1)=-dum(2,1)
188 0 : else if (nspinor==2) then ! TODO rotate wavefunction in spinor space.
189 0 : swp(:)=dum(:,1)
190 0 : dum(1,1)= dum(1,2)
191 0 : dum(2,1)=-dum(2,2)
192 0 : dum(1,2)=-swp(1)
193 0 : dum(2,2)= swp(2)
194 : end if
195 : end if
196 : !
197 : ! ==== Save values ====
198 5106882 : do ispinor=1,nspinor
199 1702294 : ibsp_bz=ibsp_bz+1
200 1702294 : Cprj_bz(iatom,ibsp_bz)%cp(1,jlmn)=dum(1,ispinor)
201 3404588 : Cprj_bz(iatom,ibsp_bz)%cp(2,jlmn)=dum(2,ispinor)
202 : end do
203 : end do !iband
204 :
205 1914462 : jl0=jl; jln0=jln
206 : end do !jlmn
207 : end do !iatom
208 :
209 104488 : call pawcprj_free(Cprjnk_kibz)
210 104488 : ABI_FREE(DS_mmpl)
211 :
212 127533 : end subroutine paw_symcprj
213 : !!***
214 :
215 : !----------------------------------------------------------------------
216 :
217 : !!****f* m_paw_sym/paw_symcprj_op
218 : !! NAME
219 : !! paw_symcprj_op
220 : !!
221 : !! FUNCTION
222 : !! Symetrize the projections cprj=<n,k|p_i> (p_i=NL PAW projector) - in-place version
223 : !!
224 : !! INPUTS
225 : !! ik_ibz=The index of the k-point in the full BZ where the matrix elements have to be symmetrized.
226 : !! nspinor=Number of spinorial components
227 : !! nband_k=Number of bands stored in cprjnk_kibz for this k-point.
228 : !! Cryst<crystal_t>=data type gathering information on unit cell and symmetries.
229 : !! %ntypat=number of type of atoms
230 : !! %natom=number of atoms in the unit cell
231 : !! %typat(natom)=type of each atom
232 : !! %indsym(4,nsym,natom)=indirect indexing array:
233 : !! for each isym,iatom, fourth element is label of atom into which iatom is sent by the INVERSE of the
234 : !! symmetry operation symrel(isym); first three elements are the primitive translations that must be subtracted
235 : !! after the transformation to get back to the original unit cell.
236 : !! Kmesh<kmesh_t>: datatype gathering information on the k-point sampling.
237 : !! %nbz=number of k-points in the full Brillouin zone
238 : !! %nibz=number of k-points in the irreducible wedge
239 : !! %tab(nkbz)=table giving for each k-point in the BZ (array kbz), the corresponding irred. point in the IBZ.
240 : !! i.e k_BZ = (IS) kIBZ where S is one of the symrec operations and I is the inversion or the identity
241 : !! %tabi(nkbz)=for each k-point in the BZ defines whether inversion has to be considered in the
242 : !! relation k_BZ=(IS) k_IBZ (1 => only S; -1 => -S)
243 : !! %tabo(nkbz)= the symmetry operation S that takes k_IBZ to each k_BZ
244 : !! Pawtab(Cryst%ntypat) <type(pawtab_type)>=paw tabulated starting data.
245 : !! Pawang <type(pawang_type)>=paw angular mesh and related data
246 : !! %lmax=Max angular momentum included in the PAW datasets used. mentioned at the second line of the psp file
247 : !! %zarot(2*lmax+1,2*lmax+1,lmax+1,nsym)=coefficients of the transformation of real spherical
248 : !! harmonics under the symmetry operations.
249 : !! in_Cprj(Cryst%natom,nspinor*nband_k)<pawcprj_type>=Input cprj
250 : !!
251 : !! OUTPUT
252 : !! out_Cprj(Cryst%natom,nspinor*nband_k)<pawcprj_type>=Symmetrized cprj matrix elements.
253 : !!
254 : !! NOTES
255 : !! Derivatives are not symmetrized.
256 : !!
257 : !! SOURCE
258 :
259 216184 : subroutine paw_symcprj_op(ik_bz,nspinor,nband_k,Cryst,Kmesh,Pawtab,Pawang,in_Cprj,out_Cprj)
260 :
261 : !Arguments ------------------------------------
262 : !scalars
263 : integer,intent(in) :: nspinor,nband_k,ik_bz
264 : type(crystal_t),intent(in) :: Cryst
265 : type(kmesh_t),intent(in) :: Kmesh
266 : type(Pawang_type),intent(in) :: Pawang
267 : !arrays
268 : type(Pawtab_type),target,intent(in) :: Pawtab(Cryst%ntypat)
269 : type(pawcprj_type),intent(in) :: in_Cprj(Cryst%natom,nspinor*nband_k)
270 : type(pawcprj_type),intent(inout) :: out_Cprj(Cryst%natom,nspinor*nband_k) !vz_i
271 :
272 : !Local variables-------------------------------
273 : !scalars
274 : integer :: iatom,iat_sym,iband,ibsp_bz
275 : integer :: ibsp_ibz,ik_ibz,indexj,ispinor,isym,itim
276 : integer :: itypat,jl,jl0,jlmn,jln,jln0,jlpm,jm,jn,lmax,mm,ncpgr
277 : real(dp) :: arg,wtk
278 : logical :: isirred
279 : !arrays
280 : integer :: r0(3) !,nlmn_atom(Cryst%natom)
281 216184 : integer,ABI_CONTIGUOUS pointer :: indlmn(:,:)
282 12572 : real(dp) :: dum(2,nspinor),kbz(3),kirr(3),phase(2),swp(2),tmp(2,nspinor)
283 216184 : real(dp),allocatable :: DS_mmpl(:,:,:)
284 :
285 : ! *********************************************************************
286 :
287 216184 : ncpgr = in_Cprj(1,1)%ncpgr
288 0 : ABI_CHECK(ncpgr==0,"Derivatives of cprj are not coded")
289 :
290 : !Get the index of the IBZ image associated to the BZ k-point ik_bz and related simmetry.
291 216184 : call Kmesh%get_BZ_item(ik_bz,kbz,ik_ibz,isym,itim,isirred=isirred)
292 :
293 216184 : if (isirred) then ! It is a point in the IBZ, Symmetrization is not needed.
294 12572 : call pawcprj_copy(in_Cprj,out_Cprj)
295 : RETURN
296 : end if
297 : !
298 : !The corresponding point kirr in the IBZ.
299 203612 : call Kmesh%get_IBZ_item(ik_ibz,kirr,wtk)
300 : !
301 : !=== DS_mmpl is the rotation matrix for real spherical harmonics associated to symrec(:,:,isym) ===
302 : !* Note the convention used by Blanco in Eq. 27 : DS_mmp multiply spherical harmonics as row vectors
303 203612 : lmax=Pawang%l_max-1 ! l_max is Max l+1
304 1221672 : ABI_MALLOC(DS_mmpl,(2*lmax+1,2*lmax+1,lmax+1))
305 5701136 : DS_mmpl=Pawang%zarot(:,:,:,isym)
306 :
307 : !Local copy.
308 : !do iatom=1,Cryst%natom
309 : !nlmn_atom(iatom)=Pawtab(Cryst%typat(iatom))%lmn_size
310 : !end do
311 : !call pawcprj_alloc(out_Cprj,ncpgr,nlmn_atom)
312 : !
313 : !===========================================
314 : !==== Loop over atoms to be symmetrized ====
315 : !===========================================
316 610836 : do iatom=1,Cryst%natom
317 407224 : itypat=Cryst%typat(iatom)
318 407224 : iat_sym=Cryst%indsym(4,isym,iatom)
319 407224 : indlmn => Pawtab(itypat)%indlmn
320 1628896 : r0=Cryst%indsym(1:3,isym,iatom) ! R^{-1} (xred(:,iatom)-tnons) = xred(:,iat_sym) + r0.
321 1628896 : arg=two_pi*dot_product(kirr,r0)
322 407224 : phase(1)=COS(arg)
323 407224 : phase(2)=SIN(arg)
324 : !
325 : ! Loop over the (jl,jm,jn) components to be symmetrized.
326 407224 : jl0=-1; jln0=-1; indexj=1
327 3863168 : do jlmn=1,Pawtab(itypat)%lmn_size
328 3252332 : jl =indlmn(1,jlmn)
329 3252332 : jm =indlmn(2,jlmn)
330 3252332 : jn =indlmn(3,jlmn)
331 3252332 : jln =indlmn(5,jlmn)
332 3252332 : jlpm=1+jl+jm
333 3252332 : if (jln/=jln0) indexj=indexj+2*jl0+1
334 : !
335 : ! === For each band, calculate contribution due to rotated real spherical harmonics ===
336 : ! FIXME check this expression; according to Blanco I should have D(S^-1} but it seems D(S) is correct
337 : ! Recheck spinorial case, presently is wrong
338 3252332 : ibsp_ibz=0
339 3252332 : ibsp_bz=0
340 6504664 : do iband=1,nband_k
341 :
342 13009328 : tmp(:,:)=zero
343 6504664 : do ispinor=1,nspinor
344 3252332 : ibsp_ibz=ibsp_ibz+1
345 14632764 : do mm=1,2*jl+1
346 8128100 : tmp(1,ispinor)=tmp(1,ispinor)+DS_mmpl(mm,jlpm,jl+1)*in_Cprj(iat_sym,ibsp_ibz)%cp(1,indexj+mm)
347 11380432 : tmp(2,ispinor)=tmp(2,ispinor)+DS_mmpl(mm,jlpm,jl+1)*in_Cprj(iat_sym,ibsp_ibz)%cp(2,indexj+mm)
348 : end do
349 : end do !ispinor
350 : !
351 : ! * Apply the phase to account if the symmetric atom belongs to a different unit cell.
352 6504664 : do ispinor=1,nspinor
353 3252332 : dum(1,ispinor)=tmp(1,ispinor)*phase(1)-tmp(2,ispinor)*phase(2)
354 6504664 : dum(2,ispinor)=tmp(1,ispinor)*phase(2)+tmp(2,ispinor)*phase(1)
355 : end do
356 : !
357 : ! * If required, apply time-reversal symmetry to retrieve the correct point in the BZ.
358 3252332 : if (itim==2) then
359 0 : if (nspinor==1) then
360 0 : dum(2,1)=-dum(2,1)
361 0 : else if (nspinor==2) then ! TODO rotate wavefunction in spinor space.
362 0 : swp(:)=dum(:,1)
363 0 : dum(1,1)= dum(1,2)
364 0 : dum(2,1)=-dum(2,2)
365 0 : dum(1,2)=-swp(1)
366 0 : dum(2,2)= swp(2)
367 : end if
368 : end if
369 : !
370 : ! ==== Save values ====
371 9756996 : do ispinor=1,nspinor
372 3252332 : ibsp_bz=ibsp_bz+1
373 3252332 : out_Cprj(iatom,ibsp_bz)%cp(1,jlmn)=dum(1,ispinor)
374 6504664 : out_Cprj(iatom,ibsp_bz)%cp(2,jlmn)=dum(2,ispinor)
375 : end do
376 : end do !iband
377 :
378 3659556 : jl0=jl; jln0=jln
379 : end do !jlmn
380 : end do !iatom
381 :
382 203612 : ABI_FREE(DS_mmpl)
383 :
384 216184 : end subroutine paw_symcprj_op
385 : !!***
386 :
387 : !----------------------------------------------------------------------
388 :
389 : END MODULE m_paw_sym
390 : !!***
|
