Line data Source code
1 : !!****m* ABINIT/m_double_grid
2 : !! NAME
3 : !! m_double_grid
4 : !!
5 : !! FUNCTION
6 : !! This module defines the double grid object. This object contains the coarse mesh
7 : !! and the dense mesh used for the interpolation of the BSE Hamiltonian,
8 : !! and contains the mapping between the two meshes.
9 : !!
10 : !! COPYRIGHT
11 : !! Copyright (C) 2008-2026 ABINIT group (YG, SP, MJV)
12 : !! This file is distributed under the terms of the
13 : !! GNU General Public License, see ~abinit/COPYING
14 : !! or http://www.gnu.org/copyleft/gpl.txt .
15 : !!
16 : !! SOURCE
17 :
18 : #if defined HAVE_CONFIG_H
19 : #include "config.h"
20 : #endif
21 :
22 : #include "abi_common.h"
23 :
24 : MODULE m_double_grid
25 :
26 : use defs_basis
27 : use m_errors
28 : use m_abicore
29 : use m_hide_blas
30 : use m_bz_mesh
31 : use m_krank
32 :
33 : use m_numeric_tools, only : wrap2_zero_one, interpol3d_indices
34 : use m_matrix, only : matr3inv
35 :
36 : implicit none
37 :
38 : private
39 : !!***
40 :
41 : !!****t* m_double_grid/double_grid_t
42 : !! NAME
43 : !! double_grid_t
44 : !!
45 : !! FUNCTION
46 : !! The double grid contains a coarse mesh and a dense mesh
47 : !! It also contains the mapping between the two meshes
48 : !!
49 : !! SOURCE
50 :
51 : type,public :: double_grid_t
52 :
53 : integer :: kmult(3)
54 : ! Number of subdivisions in the coarse box in each direction
55 :
56 : integer :: ndiv
57 : ! Total number of small sub-boxes in the large box associated to the coarse k-mesh.
58 : ! i.e. product(kmult)
59 :
60 : integer :: maxcomp_coarse(3)
61 : ! Dimensions of the box containing the coarse points in integer coord
62 :
63 : integer :: nbz_coarse
64 : ! Number of k-points in the coarse BZ (open mesh)
65 :
66 : integer :: nbz_closedcoarse
67 : ! Number of k-points inside the coarse BZ (closed mesh)
68 : ! = PROD(maxcomp_coarse+1)
69 :
70 : integer :: nbz_dense
71 : ! Number of k-point inside the dense BZ (open mesh)
72 : ! = PROD(maxcomp_coarse.*kmult)
73 :
74 : integer, allocatable :: inttoik_coarse(:)
75 : ! inttoik_coarse(nbz_closedcoarse)
76 : ! Index of the kpoint in the coarse BZ.
77 :
78 : integer, allocatable :: iktoint_coarse(:)
79 : ! iktoint_coarse(nbz_coarse)
80 : ! Index of int by kpoint
81 :
82 : integer, allocatable :: inttoik_dense(:)
83 :
84 : integer, allocatable :: iktoint_dense(:)
85 :
86 : integer,allocatable :: indices_coarse(:,:)
87 : ! indices_coarse(3,nbz_closedcoarse)
88 : ! Indices (i1,i2,i3) for each point, ordinated by the integer coord
89 :
90 : integer,allocatable :: indices_dense(:,:)
91 : ! indices_dense(6,nbz_dense)
92 : ! Indices (i1,i2,i3);(j1,j2,j3) for each point, ordinated by integer coord
93 :
94 : integer :: kptrlatt_dense(3,3)
95 : ! kptrlatt of the dense mesh
96 :
97 : real(dp) :: klatt_dense(3,3)
98 :
99 : integer :: nshiftk_dense
100 : ! Number of shifts in the dense mesh.
101 :
102 : real(dp),allocatable :: shiftk_dense(:,:)
103 : ! shift_dense(3,nshiftk_dense)
104 : ! Shifts of the dense mesh.
105 :
106 : ! Dense lattice
107 : integer :: kptrlatt_coarse(3,3)
108 : ! kptrlatt of the coarse mesh
109 :
110 : real(dp) :: klatt_coarse(3,3)
111 :
112 : integer :: nshiftk_coarse
113 : ! Number of shifts in the coarse mesh.
114 :
115 : real(dp),allocatable :: shiftk_coarse(:,:)
116 : ! shift_coarse(3,nshiftk_coarse)
117 : ! Shifts of the coarse mesh.
118 :
119 : ! Coarse lattice
120 : integer, allocatable :: g0_coarse(:,:)
121 : integer, allocatable :: g0_dense(:,:)
122 : ! g0_dense/coarse(3,nkpt_closedcoarse/dense)
123 : ! G0 vector between the kpt obtained with indices
124 : ! and the kpt obtained insize bz
125 :
126 : integer, allocatable :: dense_to_coarse(:)
127 : ! dense_to_coarse(nbz_dense)
128 : ! Give the ibz_coarse corresponding to the dense mesh (the (0,0,0) point)
129 :
130 : integer, allocatable :: coarse_to_dense(:,:)
131 : ! coarse_to_dense(nbz_coarse,ndiv)
132 : ! Give all the ibz_dense corresponding to the (0,0,0) coarse point
133 :
134 : end type double_grid_t
135 :
136 : public :: double_grid_init ! Initializes the double grid with coarse mesh and dense mesh read from file
137 : public :: double_grid_free ! Deallocate all memory
138 : public :: get_kpt_from_indices_coarse ! Returns the k-point index and g0 vector associated to the set of indices
139 : public :: compute_corresp ! Compute correspondance data between k-dense and k-coarse
140 : !public :: get_kpt_from_indices_dense
141 :
142 : public :: kptfine_av ! Find the k-points of a fine grid that are around a k-point of a coarse mesh.
143 : public :: k_neighbors ! Find 8 neighbors of given k-point on a coarse grid, and return
144 : !!***
145 :
146 : !----------------------------------------------------------------------
147 :
148 : CONTAINS !=============================================================================
149 : !!***
150 :
151 : !!****f* m_double_grid/double_grid_init
152 : !! NAME
153 : !! double_grid_init
154 : !!
155 : !! FUNCTION
156 : !! Initialize the double_grid datatype "grid" from coarse and dense mesh
157 : !!
158 : !! INPUTS
159 : !! Kmesh_coarse = descriptor of the coarse BZ sampling
160 : !! Kmesh_dense = descriptor of the dense BZ sampling
161 : !! kptrlatt_coarse(3,3) = vectors in R space that defines the reciprocal cell
162 : !! kmult(3) = multiplication factors from coarse to dense
163 : !!
164 : !! OUTPUT
165 : !! grid = double_grid to be created
166 : !!
167 : !! SOURCE
168 :
169 4 : subroutine double_grid_init(Kmesh_coarse, Kmesh_dense, kptrlatt_coarse, kmult, grid)
170 :
171 : !Argument ------------------------------------
172 : !scalars
173 : type(double_grid_t),intent(out) :: grid
174 : type(kmesh_t),intent(in) :: Kmesh_coarse,Kmesh_dense
175 : !arrays
176 : integer,intent(in) :: kptrlatt_coarse(3,3),kmult(3)
177 :
178 : !Local variables -----------------------------
179 : !scalars
180 : integer :: ii, info
181 : !arrays
182 : integer :: ipiv(3)
183 : real(dp) :: rlatt_coarse(3,3),klatt_coarse(3,3),curmat(3,3)
184 : !*********************************************
185 :
186 4 : ABI_CHECK(Kmesh_coarse%nshift == 1, "Coarse mesh works only with nshiftk=1")
187 4 : ABI_CHECK(Kmesh_dense%nshift == 1, "Dense mesh : Works only with nshiftk=1")
188 :
189 4 : grid%nshiftk_coarse = Kmesh_coarse%nshift
190 4 : grid%nshiftk_dense = Kmesh_dense%nshift
191 :
192 4 : grid%nbz_coarse = Kmesh_coarse%nbz
193 :
194 12 : ABI_MALLOC(grid%shiftk_coarse,(3,grid%nshiftk_coarse))
195 12 : ABI_MALLOC(grid%shiftk_dense,(3,grid%nshiftk_dense))
196 :
197 20 : grid%shiftk_coarse(:,:) = Kmesh_coarse%shift(:,:)
198 20 : grid%shiftk_dense(:,:) = Kmesh_dense%shift(:,:)
199 :
200 52 : grid%kptrlatt_coarse(:,:) = kptrlatt_coarse(:,:)
201 52 : rlatt_coarse(:,:) = kptrlatt_coarse(:,:)
202 4 : call matr3inv(rlatt_coarse,klatt_coarse)
203 52 : grid%klatt_coarse(:,:) = klatt_coarse(:,:)
204 :
205 4 : grid%nbz_dense = Kmesh_dense%nbz
206 4 : grid%nbz_coarse = Kmesh_coarse%nbz
207 :
208 16 : grid%kmult(:) = kmult(:)
209 4 : grid%ndiv = kmult(1)*kmult(2)*kmult(3)
210 :
211 12 : ABI_MALLOC(grid%indices_dense,(6,Kmesh_dense%nbz))
212 12 : ABI_MALLOC(grid%g0_dense,(3,Kmesh_dense%nbz))
213 12 : ABI_MALLOC(grid%iktoint_dense,(Kmesh_dense%nbz))
214 8 : ABI_MALLOC(grid%inttoik_dense,(Kmesh_dense%nbz))
215 :
216 16 : grid%maxcomp_coarse(:) = -1
217 :
218 52 : curmat(:,:) = grid%kptrlatt_coarse(:,:)
219 :
220 : ! Gaussian elimination
221 4 : call dgetrf(3,3,curmat,3,ipiv,info)
222 :
223 4 : grid%nbz_closedcoarse = 1
224 :
225 16 : do ii = 1,3
226 12 : grid%maxcomp_coarse(ii) = ABS(NINT(curmat(ipiv(ii),ipiv(ii))))
227 16 : grid%nbz_closedcoarse = grid%nbz_closedcoarse*(grid%maxcomp_coarse(ii)+1)
228 : end do
229 :
230 12 : ABI_MALLOC(grid%indices_coarse,(3,grid%nbz_closedcoarse))
231 8 : ABI_MALLOC(grid%g0_coarse,(3,grid%nbz_closedcoarse))
232 12 : ABI_MALLOC(grid%iktoint_coarse,(Kmesh_coarse%nbz))
233 12 : ABI_MALLOC(grid%inttoik_coarse,(grid%nbz_closedcoarse))
234 :
235 : ! We should pass 'grid' at this stage !
236 :
237 : call create_indices_coarse(Kmesh_coarse%bz, Kmesh_coarse%nbz, grid%klatt_coarse, &
238 : & grid%nshiftk_coarse, grid%shiftk_coarse, grid%maxcomp_coarse, grid%nbz_closedcoarse, grid%indices_coarse, &
239 4 : & grid%g0_coarse,grid%iktoint_coarse,grid%inttoik_coarse)
240 :
241 : call create_indices_dense(grid%klatt_coarse, grid%maxcomp_coarse, Kmesh_dense%bz, Kmesh_dense%nbz, &
242 : & grid%nshiftk_dense, grid%shiftk_dense, grid%kmult, grid%indices_dense, grid%g0_dense, grid%iktoint_dense, &
243 4 : & grid%inttoik_dense)
244 :
245 8 : ABI_MALLOC(grid%dense_to_coarse,(Kmesh_dense%nbz))
246 16 : ABI_MALLOC(grid%coarse_to_dense,(Kmesh_coarse%nbz,grid%ndiv))
247 :
248 : call compute_neighbours(grid%nbz_dense, grid%iktoint_dense, grid%indices_dense, &
249 : & grid%maxcomp_coarse, grid%inttoik_coarse, grid%g0_coarse, grid%nbz_closedcoarse, grid%nbz_coarse,&
250 4 : & grid%ndiv, grid%dense_to_coarse, grid%coarse_to_dense)
251 :
252 4 : end subroutine double_grid_init
253 : !!***
254 :
255 : !----------------------------------------------------------------------------
256 :
257 : !!****f* m_double_grid/create_indices_coarse
258 : !! NAME
259 : !! create_indices_coarse
260 : !!
261 : !! FUNCTION
262 : !! Create mapping between kpoints and integer indexing
263 : !!
264 : !! INPUTS
265 : !! bz(3,nbz) = k-points in the Brillouin Zone
266 : !! nbz = number of k-points
267 : !! klatt(3,3) = reciprocal space vectors defining the reciprocal cell
268 : !! nshiftk = Number of shifts
269 : !! shiftk(3,nshiftk) = Shiftks of the Brillouin Zone
270 : !! maxcomp(3) = Maximum int along each direction
271 : !! nbz_closed = Number of k-points inside the closed Brillouin Zone (adding periodic images)
272 : !!
273 : !! OUTPUT
274 : !! indices(3,nbz_closed) = indices for each k-point in the closed BZ
275 : !! g0(3,nbz_closed) = g vectors between k-point inside bz and k-point given by indices
276 : !! iktoint(nbz) = mapping between k-points in the bz and int indices
277 : !! inttoik(nbz_closed) = mapping between int indices and k-points in the bz
278 : !!
279 : !! SOURCE
280 :
281 4 : subroutine create_indices_coarse(bz, nbz, klatt, nshiftk, shiftk, maxcomp, nbz_closed, indices, g0, iktoint, inttoik)
282 :
283 : !Argument ------------------------------------
284 : !scalars
285 : integer,intent(in) :: nbz,nshiftk,nbz_closed
286 : !arrays
287 : integer,intent(in) :: maxcomp(3)
288 : integer,intent(out) :: indices(3,nbz_closed)
289 : integer,intent(out) :: g0(3,nbz_closed)
290 : integer,intent(out) :: iktoint(nbz), inttoik(nbz_closed)
291 : real(dp),intent(in) :: bz(3,nbz),klatt(3,3),shiftk(3,nshiftk)
292 :
293 : !Local variables -----------------------------
294 : !scalars
295 : integer :: ik,ii,i1,i2, i3
296 : logical :: found
297 : !arrays
298 : integer :: curg0(3)
299 : real(dp) :: curk1(3),ktoget(3)
300 : !*********************************************
301 :
302 4 : ABI_CHECK(nshiftk==1,"nshiftk != 1 not supported")
303 :
304 16 : do i1 = 0,maxcomp(1)
305 52 : do i2 = 0,maxcomp(2)
306 156 : do i3 = 0,maxcomp(3)
307 108 : ii = (i1*(maxcomp(2)+1)+i2)*(maxcomp(3)+1)+i3+1
308 432 : ktoget(:) = shiftk(:,1)+(/i1,i2,i3/)
309 1404 : curk1(:) = MATMUL(klatt(:,:),ktoget(:))
310 360 : found = .FALSE.
311 360 : do ik = 1,nbz
312 360 : if(isamek(curk1(:),bz(:,ik),curg0)) then
313 432 : indices(:,ii) = (/i1,i2,i3/)
314 432 : g0(:,ii) = curg0
315 108 : if (i1 /= maxcomp(1) .and. i2 /= maxcomp(2) .and. i3 /= maxcomp(3)) then
316 32 : iktoint(ik) = ii
317 : end if
318 108 : inttoik(ii) = ik
319 108 : found = .TRUE.
320 : exit
321 : end if
322 : end do
323 36 : if (.not. found) then
324 0 : write(std_out,*) "curk1 = ",curk1
325 0 : write(std_out,*) bz
326 0 : ABI_ERROR("A k-point generated from kptrlatt cannot be found in the BZ")
327 : end if
328 : end do
329 : end do
330 : end do
331 :
332 4 : end subroutine create_indices_coarse
333 : !!***
334 :
335 : !----------------------------------------------------------------------
336 :
337 : !!****f* m_double_grid/get_kpt_from_indices_coarse
338 : !! NAME
339 : !! get_kpt_from_indices_coarse
340 : !!
341 : !! FUNCTION
342 : !! Returns the k-point index and g0 vector associated to the set of indices
343 : !!
344 : !! INPUTS
345 : !! indices(3) = index of the searched k-point
346 : !! maxcomp(3) = Maximum int along each direction
347 : !! inttoik(nkpt) = mapping between int indices and k-points in the bz
348 : !! allg0(3,nkpt) = g vectors between k-point inside bz and k-point given by indices
349 : !! nkpt = number of k-points
350 : !!
351 : !! OUTPUT
352 : !! ikpt = index of k-point we search
353 : !! g0(3) = g-vector obtained
354 : !!
355 : !! SOURCE
356 :
357 21056 : subroutine get_kpt_from_indices_coarse(indices,maxcomp,inttoik,allg0,nkpt,ikpt,g0)
358 :
359 : !Argument ------------------------------------
360 : !scalars
361 : integer,intent(in) :: nkpt
362 : integer,intent(out) :: ikpt
363 : !arrays
364 : integer,intent(in) :: indices(3),maxcomp(3)
365 : integer,intent(in) :: inttoik(nkpt),allg0(3,nkpt)
366 : integer,intent(out) :: g0(3)
367 :
368 : !Local variables -----------------------------
369 : !scalars
370 : integer :: curicoord
371 : !*********************************************
372 :
373 21056 : curicoord = (indices(1)*(maxcomp(2)+1)+indices(2))*(maxcomp(3)+1)+indices(3)+1
374 21056 : ikpt = inttoik(curicoord)
375 84224 : g0 = allg0(:,curicoord)
376 :
377 21056 : end subroutine get_kpt_from_indices_coarse
378 : !!***
379 :
380 : !----------------------------------------------------------------------
381 :
382 : !!****f* m_double_grid/create_indices_dense
383 : !! NAME
384 : !! create_indices_dense
385 : !!
386 : !! FUNCTION
387 : !! Create mapping between kpoints and integer indexing
388 : !!
389 : !! INPUTS
390 : !! klatt_coarse(3,3) = reciprocal space vectors defining the reciprocal cell of coarse BZ
391 : !! maxcomp(3) = Maximum int along each direction
392 : !! bz_dense(3,nbz_dense) = k-points in the dense BZ
393 : !! nbz_dense = number of k-points in the dense BZ
394 : !! nshiftk = Number of shifts
395 : !! shiftk(3,nshiftk) = Shiftks of the Brillouin Zone
396 : !! kmult(3) = multiplication factors
397 : !! nbz_coarse = number of k-points in the coarse BZ
398 : !! kptrlatt_coarse(3,3) = real space vectors defining the reciprocal cell of coarse BZ
399 : !!
400 : !! OUTPUT
401 : !! indices(6,nbz_dense) = indices for each k-point in the closed BZ
402 : !! g0(3,nbz_dense) = g vectors between k-point inside bz and k-point given by indices
403 : !! iktoint(nbz_dense) = mapping between k-points in the bz and int indices
404 : !! inttoik(nbz_dense) = mapping between int indices and k-points in the bz
405 : !!
406 : !! SOURCE
407 :
408 4 : subroutine create_indices_dense(klatt_coarse, maxcomp, &
409 4 : & bz_dense, nbz_dense, nshiftk, shiftk, kmult, indices, g0, inttoik, iktoint)
410 :
411 : !Argument ------------------------------------
412 : !scalars
413 : integer,intent(in) :: nbz_dense, nshiftk
414 : !arrays
415 : integer,intent(in) :: kmult(3),maxcomp(3)
416 : integer,intent(out) :: indices(6,nbz_dense),g0(3,nbz_dense)
417 : integer,intent(out) :: inttoik(nbz_dense),iktoint(nbz_dense)
418 : real(dp),intent(in) :: bz_dense(3,nbz_dense),klatt_coarse(3,3)
419 : real(dp),intent(in) :: shiftk(3,nshiftk)
420 :
421 : !Local variables -----------------------------
422 : integer :: ik,ii,ii_coarse
423 : integer :: i1,i2,i3,j1,j2,j3
424 : logical :: found
425 : !arrays
426 : integer :: curg0(3)
427 : real(dp) :: curk1(3),ktoget(3)
428 : !*********************************************
429 :
430 4 : call wrtout(std_out, "Create Indices Dense", "COLL")
431 :
432 4 : ABI_CHECK(nshiftk==1,"nshiftk != 1 not supported")
433 :
434 12 : do i1 = 0,maxcomp(1)-1
435 28 : do i2 = 0,maxcomp(2)-1
436 56 : do i3 = 0,maxcomp(3)-1
437 32 : ii_coarse = (i1*(maxcomp(2)+1)+i2)*(maxcomp(3)+1)+i3+1
438 :
439 112 : do j1 = 0,kmult(1)-1
440 224 : do j2 = 0,kmult(2)-1
441 448 : do j3 = 0,kmult(3)-1
442 : ii = ((i1*kmult(1)+j1)*(maxcomp(2)*kmult(2)) +&
443 : & (i2*kmult(2)+j2))*(maxcomp(3)*kmult(3))+&
444 256 : & (i3*kmult(3)+j3)+1
445 :
446 256 : ktoget(1) = i1+((REAL(j1)+shiftk(1,1))/kmult(1))
447 256 : ktoget(2) = i2+((REAL(j2)+shiftk(2,1))/kmult(2))
448 256 : ktoget(3) = i3+((REAL(j3)+shiftk(3,1))/kmult(3))
449 :
450 3328 : curk1(:) = MATMUL(klatt_coarse(:,:),ktoget(:))
451 8320 : found = .FALSE.
452 8320 : do ik = 1,nbz_dense
453 8320 : if(isamek(curk1(:),bz_dense(:,ik),curg0)) then
454 1792 : indices(:,ii) = (/i1,i2,i3,j1,j2,j3/)
455 1024 : g0(:,ii) = curg0
456 256 : inttoik(ii) = ik
457 256 : iktoint(ik) = ii
458 256 : found = .TRUE.
459 : exit
460 : end if
461 : end do
462 128 : if(.not. found) then
463 0 : write(std_out,*) "curk1 = ",curk1
464 0 : write(std_out,*) bz_dense
465 0 : ABI_ERROR("Problem when creating indices")
466 : end if
467 : end do
468 : end do
469 : end do
470 :
471 : end do
472 : end do
473 : end do
474 :
475 4 : end subroutine create_indices_dense
476 : !!***
477 :
478 : !----------------------------------------------------------------------
479 :
480 : !!****f* m_double_grid/get_kpt_from_indices_dense
481 : !! NAME
482 : !! get_kpt_from_indices_coarse
483 : !!
484 : !! FUNCTION
485 : !! Returns the k-point index and g0 vector associated to the set of indices
486 : !!
487 : !! INPUTS
488 : !! indices(6) = index of the searched k-point
489 : !! maxcomp(3) = Maximum int along each direction
490 : !! kmult(3) = multiplication factors
491 : !! inttoik(nkpt) = mapping between int indices and k-points in the bz
492 : !! allg0(3,nkpt) = g vectors between k-point inside bz and k-point given by indices
493 : !! nkpt = number of k-points
494 : !!
495 : !! OUTPUT
496 : !! ikpt = index of k-point we search
497 : !! g0(3) = g-vector obtained
498 : !!
499 : !! SOURCE
500 :
501 : subroutine get_kpt_from_indices_dense(indices,maxcomp,kmult,inttoik,allg0,nkpt,ikpt,g0)
502 :
503 : !Argument ------------------------------------
504 : !scalars
505 : integer, intent(in) :: nkpt
506 : integer, intent(out) :: ikpt
507 : !arrays
508 : integer, intent(in) :: indices(6),maxcomp(3),inttoik(nkpt)
509 : integer, intent(in) :: allg0(3,nkpt),kmult(3)
510 : integer, intent(out) :: g0(3)
511 :
512 : !Local variables -----------------------------
513 : !scalars
514 : integer :: curicoord
515 : !*********************************************
516 :
517 : curicoord = ((indices(1)*kmult(1)+indices(4))*(maxcomp(2)*kmult(2))+&
518 : (indices(2)*kmult(2)+indices(5)))*(maxcomp(3)*kmult(3))+&
519 : (indices(3)*kmult(3)+indices(6))+1
520 :
521 : ikpt = inttoik(curicoord)
522 : g0 = allg0(:,curicoord)
523 :
524 : end subroutine get_kpt_from_indices_dense
525 : !!***
526 :
527 : !----------------------------------------------------------------------
528 :
529 : !!****f* m_double_grid/compute_neighbours
530 : !! NAME
531 : !! compute_neighbours
532 : !!
533 : !! FUNCTION
534 : !! Compute correspondance between points in the dense BZ and in the coarse BZ
535 : !!
536 : !! INPUTS
537 : !! nbz_dense, nbz_closedcoarse, nbz_coarse, ndiv
538 : !! iktoint_dense(nbz_dense)
539 : !! indices_dense(6,nbz_dense)
540 : !! maxcomp_coarse(3)
541 : !! inttoik_coarse(nbz_closedcoarse)
542 : !! g0_coarse(3,nbz_closedcoarse)
543 : !!
544 : !! OUTPUT
545 : !! dense_to_coarse(nbz_dense)
546 : !! coarse_to_dense(nbz_coarse,ndiv)
547 : !!
548 : !! SOURCE
549 :
550 4 : subroutine compute_neighbours(nbz_dense, iktoint_dense, indices_dense, maxcomp_coarse, &
551 4 : & inttoik_coarse, g0_coarse, nbz_closedcoarse, nbz_coarse, ndiv, dense_to_coarse, coarse_to_dense)
552 :
553 : !Argument ------------------------------------
554 : !scalars
555 : integer,intent(in) :: nbz_dense, nbz_closedcoarse, nbz_coarse, ndiv
556 : !arrays
557 : integer,intent(in) :: iktoint_dense(nbz_dense)
558 : integer,intent(in) :: indices_dense(6,nbz_dense)
559 : integer,intent(in) :: maxcomp_coarse(3)
560 : integer,intent(in) :: inttoik_coarse(nbz_closedcoarse)
561 : integer,intent(in) :: g0_coarse(3,nbz_closedcoarse)
562 : integer,intent(out) :: dense_to_coarse(nbz_dense)
563 : integer,intent(out) :: coarse_to_dense(nbz_coarse,ndiv)
564 :
565 : !Local variables -----------------------------
566 : !scalars
567 : integer :: ik_dense, iorder, ik_coarse
568 : !arrays
569 4 : integer :: curindex(nbz_coarse)
570 : integer :: curindices_dense(6), curindices_coarse(3)
571 : integer :: g0(3)
572 : !*********************************************
573 :
574 : DBG_ENTER("COLL")
575 :
576 292 : coarse_to_dense = 1
577 260 : dense_to_coarse = 1
578 :
579 36 : curindex = 1
580 260 : do ik_dense = 1, nbz_dense
581 : ! From ik_ibz in the dense mesh -> indices_dense
582 256 : iorder = iktoint_dense(ik_dense)
583 :
584 : ! From indices_dense -> indices_coarse
585 1792 : curindices_dense = indices_dense(:,iorder)
586 1024 : curindices_coarse = curindices_dense(1:3)
587 : ! From indices_coarse -> ik_ibz in the coarse mesh
588 : call get_kpt_from_indices_coarse(curindices_coarse,maxcomp_coarse,&
589 256 : & inttoik_coarse,g0_coarse,nbz_closedcoarse,ik_coarse,g0)
590 :
591 256 : dense_to_coarse(ik_dense) = ik_coarse
592 256 : coarse_to_dense(ik_coarse, curindex(ik_coarse)) = ik_dense
593 :
594 260 : curindex(ik_coarse) = curindex(ik_coarse) + 1
595 : end do
596 :
597 : DBG_EXIT("COLL")
598 :
599 4 : end subroutine compute_neighbours
600 : !!***
601 :
602 : !---------------------------------------------------------------------
603 :
604 : !!****f* m_double_grid/compute_corresp
605 : !! NAME
606 : !! compute_corresp
607 : !!
608 : !! FUNCTION
609 : !! Pre-process tables with mapping between divisions and coarse k-points
610 : !!
611 : !! INPUTS
612 : !! double_grid
613 : !!
614 : !! OUTPUT
615 : !! div2kdense(double_grid%nbz_coarse,double_grid%ndiv)
616 : !! (k_coarse,idiv) -> k_dense
617 : !! kdense2div(double_grid%nbz_dense)
618 : !! k_dense -> idiv
619 : !!
620 : !! SOURCE
621 :
622 4 : subroutine compute_corresp(double_grid, div2kdense, kdense2div)
623 :
624 : !Argument ------------------------------------
625 : !scalars
626 : class(double_grid_t),intent(in) :: double_grid
627 : !arrays
628 : integer,intent(out) :: div2kdense(double_grid%nbz_coarse,double_grid%ndiv)
629 : integer,intent(out) :: kdense2div(double_grid%nbz_dense)
630 :
631 : !Local variables -----------------------------
632 : !scalars
633 : integer :: iorder,ik_dense,ik_coarse
634 : !arrays
635 : integer :: curindices_dense(6)
636 4 : integer,allocatable :: curindex(:)
637 : !*********************************************
638 :
639 12 : ABI_MALLOC(curindex,(double_grid%nbz_coarse))
640 36 : curindex = 1
641 :
642 260 : do ik_dense = 1,double_grid%nbz_dense
643 : ! From ik_ibz in the dense mesh -> indices_dense
644 256 : iorder = double_grid%iktoint_dense(ik_dense)
645 : !g01 = double_grid%g0_dense(:,iorder)
646 :
647 : ! From indices_dense -> indices_coarse
648 : curindices_dense = double_grid%indices_dense(:,iorder)
649 :
650 256 : ik_coarse = double_grid%dense_to_coarse(ik_dense)
651 256 : div2kdense(ik_coarse,curindex(ik_coarse)) = ik_dense
652 256 : kdense2div(ik_dense) = curindex(ik_coarse)
653 :
654 260 : curindex(ik_coarse) = curindex(ik_coarse) + 1
655 : end do
656 :
657 4 : ABI_FREE(curindex)
658 :
659 4 : end subroutine compute_corresp
660 : !!***
661 :
662 : !----------------------------------------------------------------------
663 :
664 : !!****f* m_double_grid/double_grid_free
665 : !! NAME
666 : !! double_grid_free
667 : !!
668 : !! FUNCTION
669 : !! Deallocate all dynamics entities present in a double_grid structure.
670 : !!
671 : !! SOURCE
672 :
673 4 : subroutine double_grid_free(grid)
674 :
675 : !Arguments ------------------------------------
676 : class(double_grid_t),intent(inout) :: grid
677 : ! *********************************************************************
678 :
679 : !integer
680 4 : ABI_SFREE(grid%inttoik_coarse)
681 4 : ABI_SFREE(grid%inttoik_dense)
682 4 : ABI_SFREE(grid%iktoint_coarse)
683 4 : ABI_SFREE(grid%iktoint_dense)
684 4 : ABI_SFREE(grid%indices_coarse)
685 4 : ABI_SFREE(grid%indices_dense)
686 4 : ABI_SFREE(grid%g0_coarse)
687 4 : ABI_SFREE(grid%g0_dense)
688 4 : ABI_SFREE(grid%dense_to_coarse)
689 4 : ABI_SFREE(grid%coarse_to_dense)
690 :
691 : !real
692 4 : ABI_SFREE(grid%shiftk_dense)
693 4 : ABI_SFREE(grid%shiftk_coarse)
694 :
695 4 : end subroutine double_grid_free
696 : !!***
697 :
698 : !----------------------------------------------------------------------
699 :
700 : !!****f* m_double_grid/kptfine_av
701 : !! NAME
702 : !! kptfine_av
703 : !!
704 : !! FUNCTION
705 : !! Find the k-points of a fine grid that are around a k-point of a coarse mesh.
706 : !!
707 : !! INPUTS
708 : !! center(3) = the point of the coarse mesh around which you want know which
709 : !! k-points of the fine mesh belong to.
710 : !! qptrlatt(3,3) = qptrlatt of the considered calculation (this is obtained
711 : !! from the input variable ngqpt and shiftq.
712 : !! kpt_fine(3,nkpt_fine) = this table contain all the k-points of the fine grid
713 : !! in the full BZ (no sym op. allowed) and is read from the header
714 : !! of the dense WF file.
715 : !! nkpt_fine = number of k-points of the fine grid read from the header of the
716 : !! dense WF file.
717 : !!
718 : !! OUTPUT
719 : !! kpt_fine_sub(nkpt_sub) = k-points of the fine grid that are around center(3)
720 : !! nkpt_sub = number of k-points of the fine grid that are around center(3)
721 : !! wgt_sub(nkpt_sub) = weight of the k-points of the fine grid that are around center(3).
722 : !!
723 : !! SOURCE
724 :
725 24 : subroutine kptfine_av(center,qptrlatt,kpt_fine,nkpt_fine,kpt_fine_sub,nkpt_sub,wgt_sub)
726 :
727 : !Arguments ------------------------------------
728 : !scalars
729 : integer,intent(in) :: nkpt_fine
730 : integer,intent(out) :: nkpt_sub
731 : !arrays
732 : integer,intent(in) :: qptrlatt(3,3)
733 : real(dp),intent(in) :: kpt_fine(3,nkpt_fine)
734 : real(dp),intent(in) :: center(3)
735 : integer,pointer :: kpt_fine_sub(:)
736 : real(dp),pointer :: wgt_sub(:)
737 :
738 : !Local variables-------------------------------
739 : !scalars
740 : integer :: ikpt,aa,bb,cc
741 : integer :: ii,jj
742 : !arrays
743 : real(dp) :: center_ref(3)
744 : real(dp) :: kpt_fine_ref(3)
745 : real(dp) :: kpt_tmp(3),kpt_tmp2(3)
746 24 : integer,allocatable :: kpt_fine_sub_tmp(:)
747 24 : real(dp),allocatable :: wgt_sub_tmp(:)
748 : logical :: found(3)
749 : ! *************************************************************************
750 :
751 72 : ABI_MALLOC(kpt_fine_sub_tmp,(nkpt_fine))
752 72 : ABI_MALLOC(wgt_sub_tmp,(nkpt_fine))
753 :
754 : !It is easier to work in real space using the qptrlatt matrices because in this
755 : !referential any k-points sampling will be cast into an orthorhombic shape.
756 : !In that space we can simply take all k-points of the fine grid that between
757 : !center_ref-0.5 and center_ref+0.5
758 :
759 600 : center_ref = MATMUL(qptrlatt,center)
760 :
761 : !When considering points center(3) that lying close or on a BZ edge we need to
762 : !take the k-points of the fine grid taking into account unklamp vectors. This
763 : !is done with the aa, bb and cc loops.
764 :
765 : ii = 1
766 1560 : do ikpt=1,nkpt_fine
767 6144 : kpt_tmp = kpt_fine(:,ikpt)
768 6168 : do aa=-1,1
769 4608 : kpt_tmp2(1) = kpt_tmp(1)+aa
770 19968 : do bb=-1,1
771 13824 : kpt_tmp2(2) = kpt_tmp(2)+bb
772 59904 : do cc=-1,1
773 41472 : kpt_tmp2(3) = kpt_tmp(3)+cc
774 1036800 : kpt_fine_ref = MATMUL(qptrlatt,kpt_tmp2)
775 41472 : if((kpt_fine_ref(1)>=center_ref(1)-0.5-tol8).and.&
776 13824 : & (kpt_fine_ref(1)<=center_ref(1)+0.5+tol8)) then
777 10368 : if((kpt_fine_ref(2)>=center_ref(2)-0.5-tol8).and.&
778 : & (kpt_fine_ref(2)<=center_ref(2)+0.5+tol8)) then
779 2592 : if((kpt_fine_ref(3)>=center_ref(3)-0.5-tol8).and.&
780 : & (kpt_fine_ref(3)<=center_ref(3)+0.5+tol8)) then
781 648 : kpt_fine_sub_tmp(ii) = ikpt
782 648 : ii = ii +1
783 : end if
784 : end if
785 : end if
786 : end do
787 : end do
788 : end do
789 : end do
790 :
791 24 : nkpt_sub = ii-1
792 72 : ABI_MALLOC(kpt_fine_sub,(nkpt_sub))
793 72 : ABI_MALLOC(wgt_sub,(nkpt_sub))
794 :
795 672 : do jj=1,nkpt_sub
796 672 : kpt_fine_sub(jj) = kpt_fine_sub_tmp(jj)
797 : end do
798 :
799 : !We then compute a weight function. This weight function is simply a
800 : !rectangular weight function that take the value 1 for k-points of the fine
801 : !grid inside the cube, 0.5 for k-points that are lying on one face of the cube,
802 : !0.25 for k-points that are lying on an edge of the cube and 0.125 for k-points
803 : !that are lying on a peak of the cube.
804 :
805 672 : wgt_sub(:) = 1.0
806 :
807 672 : do ikpt=1,nkpt_sub
808 2592 : found(:) = .True.
809 2592 : kpt_tmp = kpt_fine(:,kpt_fine_sub(ikpt))
810 2616 : do aa=-1,1
811 1944 : kpt_tmp2(1) = kpt_tmp(1)+aa
812 8424 : do bb=-1,1
813 5832 : kpt_tmp2(2) = kpt_tmp(2)+bb
814 25272 : do cc=-1,1
815 17496 : kpt_tmp2(3) = kpt_tmp(3)+cc
816 437400 : kpt_fine_ref = MATMUL(qptrlatt,kpt_tmp2)
817 17496 : if((ABS(kpt_fine_ref(1)-center_ref(1)-0.5)< tol8) .or.&
818 : (ABS(kpt_fine_ref(1)-center_ref(1)+0.5) < tol8)) then
819 3888 : if(found(1)) then
820 432 : wgt_sub(ikpt) = wgt_sub(ikpt)*0.5; found(1) = .False.
821 : end if
822 : end if
823 17496 : if((ABS(kpt_fine_ref(2)-center_ref(2)-0.5) < tol8) .or.&
824 : (ABS(kpt_fine_ref(2)-center_ref(2)+0.5) < tol8)) then
825 3888 : if(found(2)) then
826 432 : wgt_sub(ikpt) = wgt_sub(ikpt)*0.5; found(2) = .False.
827 : end if
828 : end if
829 17496 : if((ABS(kpt_fine_ref(3)-center_ref(3)-0.5)< tol8) .or.&
830 5832 : (ABS(kpt_fine_ref(3)-center_ref(3)+0.5) < tol8)) then
831 3888 : if(found(3)) then
832 432 : wgt_sub(ikpt) = wgt_sub(ikpt)*0.5; found(3) = .False.
833 : end if
834 : end if
835 : end do
836 : end do
837 : end do
838 : end do
839 :
840 24 : ABI_FREE(kpt_fine_sub_tmp)
841 24 : ABI_FREE(wgt_sub_tmp)
842 :
843 24 : end subroutine kptfine_av
844 : !!***
845 :
846 : !!****f* m_double_grid/k_neighbors
847 : !!
848 : !! NAME
849 : !! k_neighbors
850 : !!
851 : !! FUNCTION
852 : !! find 8 neighbors of given k-point on a coarse grid, and return
853 : !! them along with relative k-shift within coarse grid cell
854 : !!
855 : !! INPUTS
856 : !! kpt = k-point to be interpolated to, in full BZ
857 : !! kptrlatt = lattice vectors for coarse k-grid
858 : !! invrankkpt = rank list to find k-points
859 : !!
860 : !! OUTPUT
861 : !! rel_kpt = k-point coordinates renormalized to coarse grid cell
862 : !! kpt_phon_indices = indices of k-points on corners of cell
863 : !!
864 : !! TODO
865 : !! This routine is not used anymore. Deprecate or Remove?
866 : !!
867 : !! SOURCE
868 :
869 0 : subroutine k_neighbors(kpt, kptrlatt,krank, rel_kpt, kpt_phon_indices)
870 :
871 : ! inputs
872 : real(dp), intent(in) :: kpt(3)
873 : integer, intent(in) :: kptrlatt(3,3)
874 : type(krank_t), intent(in) :: krank
875 :
876 : ! outputs
877 : real(dp), intent(out) :: rel_kpt(3)
878 : integer, intent(out) :: kpt_phon_indices(8)
879 : ! local vars
880 : integer :: symrankkpt
881 : integer :: ir1,ir2,ir3, pr1,pr2,pr3
882 : real(dp) :: redkpt(3), cornerkpt(3), res
883 : ! *************************************************************************
884 :
885 : !wrap fine kpt to [0,1]
886 0 : call wrap2_zero_one(kpt(1),redkpt(1),res)
887 0 : call wrap2_zero_one(kpt(2),redkpt(2),res)
888 0 : call wrap2_zero_one(kpt(3),redkpt(3),res)
889 : !find 8 indices of points neighboring ikpt_phon, for interpolation
890 : call interpol3d_indices (redkpt,kptrlatt(1,1),kptrlatt(2,2),kptrlatt(3,3), &
891 0 : ir1,ir2,ir3, pr1,pr2,pr3)
892 :
893 : !transpose ir pr to ikpt_phon indices
894 : !order of kpt_phons:
895 : !ir1 ir2 ir3
896 0 : cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/)
897 0 : symrankkpt = krank%get_rank(cornerkpt)
898 0 : kpt_phon_indices(1) = krank%invrank(symrankkpt)
899 : !pr1 ir2 ir3
900 0 : cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/)
901 0 : symrankkpt = krank%get_rank (cornerkpt)
902 0 : kpt_phon_indices(2) = krank%invrank(symrankkpt)
903 : !ir1 pr2 ir3
904 0 : cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/)
905 0 : symrankkpt = krank%get_rank (cornerkpt)
906 0 : kpt_phon_indices(3) = krank%invrank(symrankkpt)
907 : !pr1 pr2 ir3
908 0 : cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/)
909 0 : symrankkpt = krank%get_rank (cornerkpt)
910 0 : kpt_phon_indices(4) = krank%invrank(symrankkpt)
911 : !ir1 ir2 pr3
912 0 : cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/)
913 0 : symrankkpt = krank%get_rank (cornerkpt)
914 0 : kpt_phon_indices(5) = krank%invrank(symrankkpt)
915 : !pr1 ir2 pr3
916 0 : cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/)
917 0 : symrankkpt = krank%get_rank (cornerkpt)
918 0 : kpt_phon_indices(6) = krank%invrank(symrankkpt)
919 : !ir1 pr2 pr3
920 0 : cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/)
921 0 : symrankkpt = krank%get_rank (cornerkpt)
922 0 : kpt_phon_indices(7) = krank%invrank(symrankkpt)
923 : !pr1 pr2 pr3
924 0 : cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/)
925 0 : symrankkpt = krank%get_rank (cornerkpt)
926 0 : kpt_phon_indices(8) = krank%invrank(symrankkpt)
927 :
928 : !retrieve the gkq matrix for all q, at the neighbor k vectors
929 0 : rel_kpt(1) = redkpt(1)*kptrlatt(1,1)-real(ir1-1)
930 0 : rel_kpt(2) = redkpt(2)*kptrlatt(2,2)-real(ir2-1)
931 0 : rel_kpt(3) = redkpt(3)*kptrlatt(3,3)-real(ir3-1)
932 :
933 0 : end subroutine k_neighbors
934 : !!***
935 :
936 0 : END MODULE m_double_grid
937 : !!***
|
