TABLE OF CONTENTS
- ABINIT/m_double_grid
- m_double_grid/compute_corresp
- m_double_grid/compute_neighbours
- m_double_grid/create_indices_coarse
- m_double_grid/create_indices_dense
- m_double_grid/double_grid_free
- m_double_grid/double_grid_init
- m_double_grid/double_grid_t
- m_double_grid/get_kpt_from_indices_coarse
- m_double_grid/get_kpt_from_indices_dense
- m_double_grid/k_neighbors
- m_double_grid/kptfine_av
ABINIT/m_double_grid [ Modules ]
NAME
m_double_grid
FUNCTION
This module defines the double grid object. This object contains the coarse mesh and the dense mesh used for the interpolation of the BSE Hamiltonian, and contains the mapping between the two meshes.
COPYRIGHT
Copyright (C) 2008-2024 ABINIT group (YG, SP, MJV) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
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_symtk, only : matr3inv 35 36 implicit none 37 38 private
m_double_grid/compute_corresp [ Functions ]
[ Top ] [ m_double_grid ] [ Functions ]
NAME
compute_corresp
FUNCTION
Pre-process tables with mapping between divisions and coarse k-points
INPUTS
double_grid
OUTPUT
div2kdense(double_grid%nbz_coarse,double_grid%ndiv) (k_coarse,idiv) -> k_dense kdense2div(double_grid%nbz_dense) k_dense -> idiv
SOURCE
628 subroutine compute_corresp(double_grid, div2kdense, kdense2div) 629 630 !Argument ------------------------------------ 631 !scalars 632 type(double_grid_t),intent(in) :: double_grid 633 !arrays 634 integer,intent(out) :: div2kdense(double_grid%nbz_coarse,double_grid%ndiv) 635 integer,intent(out) :: kdense2div(double_grid%nbz_dense) 636 637 !Local variables ----------------------------- 638 !scalars 639 integer :: iorder,ik_dense,ik_coarse 640 !arrays 641 integer :: curindices_dense(6) 642 integer,allocatable :: curindex(:) 643 644 !********************************************* 645 ABI_MALLOC(curindex,(double_grid%nbz_coarse)) 646 curindex = 1 647 648 do ik_dense = 1,double_grid%nbz_dense 649 650 ! From ik_ibz in the dense mesh -> indices_dense 651 iorder = double_grid%iktoint_dense(ik_dense) 652 !g01 = double_grid%g0_dense(:,iorder) 653 654 ! From indices_dense -> indices_coarse 655 curindices_dense = double_grid%indices_dense(:,iorder) 656 657 ik_coarse = double_grid%dense_to_coarse(ik_dense) 658 div2kdense(ik_coarse,curindex(ik_coarse)) = ik_dense 659 kdense2div(ik_dense) = curindex(ik_coarse) 660 661 curindex(ik_coarse) = curindex(ik_coarse) + 1 662 663 end do 664 665 ABI_FREE(curindex) 666 667 end subroutine compute_corresp
m_double_grid/compute_neighbours [ Functions ]
[ Top ] [ m_double_grid ] [ Functions ]
NAME
compute_neighbours
FUNCTION
Compute correspondance between points in the dense BZ and in the coarse BZ
INPUTS
nbz_dense, nbz_closedcoarse, nbz_coarse, ndiv iktoint_dense(nbz_dense) indices_dense(6,nbz_dense) maxcomp_coarse(3) inttoik_coarse(nbz_closedcoarse) g0_coarse(3,nbz_closedcoarse)
OUTPUT
dense_to_coarse(nbz_dense) coarse_to_dense(nbz_coarse,ndiv)
SOURCE
555 subroutine compute_neighbours(nbz_dense, iktoint_dense, indices_dense, maxcomp_coarse, & 556 & inttoik_coarse, g0_coarse, nbz_closedcoarse, nbz_coarse, ndiv, dense_to_coarse, coarse_to_dense) 557 558 !Argument ------------------------------------ 559 !scalars 560 integer,intent(in) :: nbz_dense, nbz_closedcoarse, nbz_coarse, ndiv 561 !arrays 562 integer,intent(in) :: iktoint_dense(nbz_dense) 563 integer,intent(in) :: indices_dense(6,nbz_dense) 564 integer,intent(in) :: maxcomp_coarse(3) 565 integer,intent(in) :: inttoik_coarse(nbz_closedcoarse) 566 integer,intent(in) :: g0_coarse(3,nbz_closedcoarse) 567 integer,intent(out) :: dense_to_coarse(nbz_dense) 568 integer,intent(out) :: coarse_to_dense(nbz_coarse,ndiv) 569 570 !Local variables ----------------------------- 571 !scalars 572 integer :: ik_dense, iorder, ik_coarse 573 !arrays 574 integer :: curindex(nbz_coarse) 575 integer :: curindices_dense(6), curindices_coarse(3) 576 integer :: g0(3) 577 578 !********************************************* 579 580 DBG_ENTER("COLL") 581 582 coarse_to_dense = 1 583 dense_to_coarse = 1 584 585 curindex = 1 586 do ik_dense = 1, nbz_dense 587 ! From ik_ibz in the dense mesh -> indices_dense 588 iorder = iktoint_dense(ik_dense) 589 590 ! From indices_dense -> indices_coarse 591 curindices_dense = indices_dense(:,iorder) 592 curindices_coarse = curindices_dense(1:3) 593 ! From indices_coarse -> ik_ibz in the coarse mesh 594 call get_kpt_from_indices_coarse(curindices_coarse,maxcomp_coarse,& 595 & inttoik_coarse,g0_coarse,nbz_closedcoarse,ik_coarse,g0) 596 597 dense_to_coarse(ik_dense) = ik_coarse 598 coarse_to_dense(ik_coarse, curindex(ik_coarse)) = ik_dense 599 600 curindex(ik_coarse) = curindex(ik_coarse) + 1 601 end do 602 603 DBG_EXIT("COLL") 604 605 end subroutine compute_neighbours
m_double_grid/create_indices_coarse [ Functions ]
[ Top ] [ m_double_grid ] [ Functions ]
NAME
create_indices_coarse
FUNCTION
Create mapping between kpoints and integer indexing
INPUTS
bz(3,nbz) = k-points in the Brillouin Zone nbz = number of k-points klatt(3,3) = reciprocal space vectors defining the reciprocal cell nshiftk = Number of shifts shiftk(3,nshiftk) = Shiftks of the Brillouin Zone maxcomp(3) = Maximum int along each direction nbz_closed = Number of k-points inside the closed Brillouin Zone (adding periodic images)
OUTPUT
indices(3,nbz_closed) = indices for each k-point in the closed BZ g0(3,nbz_closed) = g vectors between k-point inside bz and k-point given by indices iktoint(nbz) = mapping between k-points in the bz and int indices inttoik(nbz_closed) = mapping between int indices and k-points in the bz
SOURCE
282 subroutine create_indices_coarse(bz, nbz, klatt, nshiftk, shiftk, maxcomp, nbz_closed, indices, g0, iktoint, inttoik) 283 284 !Argument ------------------------------------ 285 !scalars 286 integer,intent(in) :: nbz,nshiftk,nbz_closed 287 !arrays 288 integer,intent(in) :: maxcomp(3) 289 integer,intent(out) :: indices(3,nbz_closed) 290 integer,intent(out) :: g0(3,nbz_closed) 291 integer,intent(out) :: iktoint(nbz), inttoik(nbz_closed) 292 real(dp),intent(in) :: bz(3,nbz),klatt(3,3),shiftk(3,nshiftk) 293 294 !Local variables ----------------------------- 295 !scalars 296 integer :: ik,ii,i1,i2, i3 297 logical :: found 298 !arrays 299 integer :: curg0(3) 300 real(dp) :: curk1(3),ktoget(3) 301 302 !********************************************* 303 304 ABI_CHECK(nshiftk==1,"nshiftk != 1 not supported") 305 306 do i1 = 0,maxcomp(1) 307 do i2 = 0,maxcomp(2) 308 do i3 = 0,maxcomp(3) 309 ii = (i1*(maxcomp(2)+1)+i2)*(maxcomp(3)+1)+i3+1 310 ktoget(:) = shiftk(:,1)+(/i1,i2,i3/) 311 curk1(:) = MATMUL(klatt(:,:),ktoget(:)) 312 found = .FALSE. 313 do ik = 1,nbz 314 if(isamek(curk1(:),bz(:,ik),curg0)) then 315 indices(:,ii) = (/i1,i2,i3/) 316 g0(:,ii) = curg0 317 if (i1 /= maxcomp(1) .and. i2 /= maxcomp(2) .and. i3 /= maxcomp(3)) then 318 iktoint(ik) = ii 319 end if 320 inttoik(ii) = ik 321 found = .TRUE. 322 exit 323 end if 324 end do 325 if (.not. found) then 326 write(std_out,*) "curk1 = ",curk1 327 write(std_out,*) bz 328 ABI_ERROR("A k-point generated from kptrlatt cannot be found in the BZ") 329 end if 330 end do 331 end do 332 end do 333 334 end subroutine create_indices_coarse
m_double_grid/create_indices_dense [ Functions ]
[ Top ] [ m_double_grid ] [ Functions ]
NAME
create_indices_dense
FUNCTION
Create mapping between kpoints and integer indexing
INPUTS
klatt_coarse(3,3) = reciprocal space vectors defining the reciprocal cell of coarse BZ maxcomp(3) = Maximum int along each direction bz_dense(3,nbz_dense) = k-points in the dense BZ nbz_dense = number of k-points in the dense BZ nshiftk = Number of shifts shiftk(3,nshiftk) = Shiftks of the Brillouin Zone kmult(3) = multiplication factors nbz_coarse = number of k-points in the coarse BZ kptrlatt_coarse(3,3) = real space vectors defining the reciprocal cell of coarse BZ
OUTPUT
indices(6,nbz_dense) = indices for each k-point in the closed BZ g0(3,nbz_dense) = g vectors between k-point inside bz and k-point given by indices iktoint(nbz_dense) = mapping between k-points in the bz and int indices inttoik(nbz_dense) = mapping between int indices and k-points in the bz
SOURCE
411 subroutine create_indices_dense(klatt_coarse, maxcomp, & 412 & bz_dense, nbz_dense, nshiftk, shiftk, kmult, indices, g0, inttoik, iktoint) 413 414 !Argument ------------------------------------ 415 !scalars 416 integer,intent(in) :: nbz_dense, nshiftk 417 !arrays 418 integer,intent(in) :: kmult(3),maxcomp(3) 419 integer,intent(out) :: indices(6,nbz_dense),g0(3,nbz_dense) 420 integer,intent(out) :: inttoik(nbz_dense),iktoint(nbz_dense) 421 real(dp),intent(in) :: bz_dense(3,nbz_dense),klatt_coarse(3,3) 422 real(dp),intent(in) :: shiftk(3,nshiftk) 423 424 !Local variables ----------------------------- 425 integer :: ik,ii,ii_coarse 426 integer :: i1,i2,i3,j1,j2,j3 427 logical :: found 428 !arrays 429 integer :: curg0(3) 430 real(dp) :: curk1(3),ktoget(3) 431 432 !********************************************* 433 434 call wrtout(std_out, "Create Indices Dense", "COLL") 435 436 ABI_CHECK(nshiftk==1,"nshiftk != 1 not supported") 437 438 do i1 = 0,maxcomp(1)-1 439 do i2 = 0,maxcomp(2)-1 440 do i3 = 0,maxcomp(3)-1 441 ii_coarse = (i1*(maxcomp(2)+1)+i2)*(maxcomp(3)+1)+i3+1 442 443 do j1 = 0,kmult(1)-1 444 do j2 = 0,kmult(2)-1 445 do j3 = 0,kmult(3)-1 446 ii = ((i1*kmult(1)+j1)*(maxcomp(2)*kmult(2)) +& 447 & (i2*kmult(2)+j2))*(maxcomp(3)*kmult(3))+& 448 & (i3*kmult(3)+j3)+1 449 450 ktoget(1) = i1+((REAL(j1)+shiftk(1,1))/kmult(1)) 451 ktoget(2) = i2+((REAL(j2)+shiftk(2,1))/kmult(2)) 452 ktoget(3) = i3+((REAL(j3)+shiftk(3,1))/kmult(3)) 453 454 curk1(:) = MATMUL(klatt_coarse(:,:),ktoget(:)) 455 found = .FALSE. 456 do ik = 1,nbz_dense 457 if(isamek(curk1(:),bz_dense(:,ik),curg0)) then 458 indices(:,ii) = (/i1,i2,i3,j1,j2,j3/) 459 g0(:,ii) = curg0 460 inttoik(ii) = ik 461 iktoint(ik) = ii 462 found = .TRUE. 463 exit 464 end if 465 end do 466 if(.not. found) then 467 write(std_out,*) "curk1 = ",curk1 468 write(std_out,*) bz_dense 469 ABI_ERROR("Problem when creating indices") 470 end if 471 end do 472 end do 473 end do 474 475 end do 476 end do 477 end do 478 479 end subroutine create_indices_dense
m_double_grid/double_grid_free [ Functions ]
[ Top ] [ m_double_grid ] [ Functions ]
NAME
double_grid_free
FUNCTION
Deallocate all dynamics entities present in a double_grid structure.
INPUTS
grid<double_grid>=The datatype to be freed.
SIDE EFFECTS
All allocated memory is released.
SOURCE
687 subroutine double_grid_free(grid) 688 689 !Arguments ------------------------------------ 690 type(double_grid_t),intent(inout) :: grid 691 692 ! ********************************************************************* 693 694 !integer 695 ABI_SFREE(grid%inttoik_coarse) 696 ABI_SFREE(grid%inttoik_dense) 697 ABI_SFREE(grid%iktoint_coarse) 698 ABI_SFREE(grid%iktoint_dense) 699 ABI_SFREE(grid%indices_coarse) 700 ABI_SFREE(grid%indices_dense) 701 ABI_SFREE(grid%g0_coarse) 702 ABI_SFREE(grid%g0_dense) 703 ABI_SFREE(grid%dense_to_coarse) 704 ABI_SFREE(grid%coarse_to_dense) 705 706 !real 707 ABI_SFREE(grid%shiftk_dense) 708 ABI_SFREE(grid%shiftk_coarse) 709 710 end subroutine double_grid_free
m_double_grid/double_grid_init [ Functions ]
[ Top ] [ m_double_grid ] [ Functions ]
NAME
double_grid_init
FUNCTION
Initialize the double_grid datatype "grid" from coarse and dense mesh
INPUTS
Kmesh_coarse = descriptor of the coarse BZ sampling Kmesh_dense = descriptor of the dense BZ sampling kptrlatt_coarse(3,3) = vectors in R space that defines the reciprocal cell kmult(3) = multiplication factors from coarse to dense
OUTPUT
grid = double_grid to be created
SOURCE
169 subroutine double_grid_init(Kmesh_coarse,Kmesh_dense,kptrlatt_coarse,kmult,grid) 170 171 !Argument ------------------------------------ 172 !scalars 173 type(kmesh_t),intent(in) :: Kmesh_coarse,Kmesh_dense 174 type(double_grid_t),intent(out) :: grid 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 187 ABI_CHECK(Kmesh_coarse%nshift == 1, "Coarse mesh works only with nshiftk=1") 188 ABI_CHECK(Kmesh_dense%nshift == 1, "Dense mesh : Works only with nshiftk=1") 189 190 grid%nshiftk_coarse = Kmesh_coarse%nshift 191 grid%nshiftk_dense = Kmesh_dense%nshift 192 193 grid%nbz_coarse = Kmesh_coarse%nbz 194 195 ABI_MALLOC(grid%shiftk_coarse,(3,grid%nshiftk_coarse)) 196 ABI_MALLOC(grid%shiftk_dense,(3,grid%nshiftk_dense)) 197 198 grid%shiftk_coarse(:,:) = Kmesh_coarse%shift(:,:) 199 grid%shiftk_dense(:,:) = Kmesh_dense%shift(:,:) 200 201 grid%kptrlatt_coarse(:,:) = kptrlatt_coarse(:,:) 202 rlatt_coarse(:,:) = kptrlatt_coarse(:,:) 203 call matr3inv(rlatt_coarse,klatt_coarse) 204 grid%klatt_coarse(:,:) = klatt_coarse(:,:) 205 206 grid%nbz_dense = Kmesh_dense%nbz 207 grid%nbz_coarse = Kmesh_coarse%nbz 208 209 grid%kmult(:) = kmult(:) 210 grid%ndiv = kmult(1)*kmult(2)*kmult(3) 211 212 ABI_MALLOC(grid%indices_dense,(6,Kmesh_dense%nbz)) 213 ABI_MALLOC(grid%g0_dense,(3,Kmesh_dense%nbz)) 214 ABI_MALLOC(grid%iktoint_dense,(Kmesh_dense%nbz)) 215 ABI_MALLOC(grid%inttoik_dense,(Kmesh_dense%nbz)) 216 217 grid%maxcomp_coarse(:) = -1 218 219 curmat(:,:) = grid%kptrlatt_coarse(:,:) 220 221 ! Gaussian elimination 222 call dgetrf(3,3,curmat,3,ipiv,info) 223 224 grid%nbz_closedcoarse = 1 225 226 do ii = 1,3 227 grid%maxcomp_coarse(ii) = ABS(NINT(curmat(ipiv(ii),ipiv(ii)))) 228 grid%nbz_closedcoarse = grid%nbz_closedcoarse*(grid%maxcomp_coarse(ii)+1) 229 end do 230 231 ABI_MALLOC(grid%indices_coarse,(3,grid%nbz_closedcoarse)) 232 ABI_MALLOC(grid%g0_coarse,(3,grid%nbz_closedcoarse)) 233 ABI_MALLOC(grid%iktoint_coarse,(Kmesh_coarse%nbz)) 234 ABI_MALLOC(grid%inttoik_coarse,(grid%nbz_closedcoarse)) 235 236 ! We should pass 'grid' at this stage ! 237 238 call create_indices_coarse(Kmesh_coarse%bz, Kmesh_coarse%nbz, grid%klatt_coarse, & 239 & grid%nshiftk_coarse, grid%shiftk_coarse, grid%maxcomp_coarse, grid%nbz_closedcoarse, grid%indices_coarse, & 240 & grid%g0_coarse,grid%iktoint_coarse,grid%inttoik_coarse) 241 242 call create_indices_dense(grid%klatt_coarse, grid%maxcomp_coarse, Kmesh_dense%bz, Kmesh_dense%nbz, & 243 & grid%nshiftk_dense, grid%shiftk_dense, grid%kmult, grid%indices_dense, grid%g0_dense, grid%iktoint_dense, & 244 & grid%inttoik_dense) 245 246 ABI_MALLOC(grid%dense_to_coarse,(Kmesh_dense%nbz)) 247 ABI_MALLOC(grid%coarse_to_dense,(Kmesh_coarse%nbz,grid%ndiv)) 248 249 call compute_neighbours(grid%nbz_dense, grid%iktoint_dense, grid%indices_dense, & 250 & grid%maxcomp_coarse, grid%inttoik_coarse, grid%g0_coarse, grid%nbz_closedcoarse, grid%nbz_coarse,& 251 & grid%ndiv, grid%dense_to_coarse, grid%coarse_to_dense) 252 253 end subroutine double_grid_init
m_double_grid/double_grid_t [ Types ]
[ Top ] [ m_double_grid ] [ Types ]
NAME
double_grid_t
FUNCTION
The double grid contains a coarse mesh and a dense mesh It also contains the mapping between the two meshes
SOURCE
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
m_double_grid/get_kpt_from_indices_coarse [ Functions ]
[ Top ] [ m_double_grid ] [ Functions ]
NAME
get_kpt_from_indices_coarse
FUNCTION
Returns the k-point index and g0 vector associated to the set of indices
INPUTS
indices(3) = index of the searched k-point maxcomp(3) = Maximum int along each direction inttoik(nkpt) = mapping between int indices and k-points in the bz allg0(3,nkpt) = g vectors between k-point inside bz and k-point given by indices nkpt = number of k-points
OUTPUT
ikpt = index of k-point we search g0(3) = g-vector obtained
SOURCE
359 subroutine get_kpt_from_indices_coarse(indices,maxcomp,inttoik,allg0,nkpt,ikpt,g0) 360 361 !Argument ------------------------------------ 362 !scalars 363 integer,intent(in) :: nkpt 364 integer,intent(out) :: ikpt 365 !arrays 366 integer,intent(in) :: indices(3),maxcomp(3) 367 integer,intent(in) :: inttoik(nkpt),allg0(3,nkpt) 368 integer,intent(out) :: g0(3) 369 370 !Local variables ----------------------------- 371 !scalars 372 integer :: curicoord 373 374 !********************************************* 375 376 curicoord = (indices(1)*(maxcomp(2)+1)+indices(2))*(maxcomp(3)+1)+indices(3)+1 377 ikpt = inttoik(curicoord) 378 g0 = allg0(:,curicoord) 379 380 end subroutine get_kpt_from_indices_coarse
m_double_grid/get_kpt_from_indices_dense [ Functions ]
[ Top ] [ m_double_grid ] [ Functions ]
NAME
get_kpt_from_indices_coarse
FUNCTION
Returns the k-point index and g0 vector associated to the set of indices
INPUTS
indices(6) = index of the searched k-point maxcomp(3) = Maximum int along each direction kmult(3) = multiplication factors inttoik(nkpt) = mapping between int indices and k-points in the bz allg0(3,nkpt) = g vectors between k-point inside bz and k-point given by indices nkpt = number of k-points
OUTPUT
ikpt = index of k-point we search g0(3) = g-vector obtained
SOURCE
505 subroutine get_kpt_from_indices_dense(indices,maxcomp,kmult,inttoik,allg0,nkpt,ikpt,g0) 506 507 !Argument ------------------------------------ 508 !scalars 509 integer, intent(in) :: nkpt 510 integer, intent(out) :: ikpt 511 !arrays 512 integer, intent(in) :: indices(6),maxcomp(3),inttoik(nkpt) 513 integer, intent(in) :: allg0(3,nkpt),kmult(3) 514 integer, intent(out) :: g0(3) 515 516 !Local variables ----------------------------- 517 !scalars 518 integer :: curicoord 519 520 !********************************************* 521 522 curicoord = ((indices(1)*kmult(1)+indices(4))*(maxcomp(2)*kmult(2))+& 523 & (indices(2)*kmult(2)+indices(5)))*(maxcomp(3)*kmult(3))+& 524 & (indices(3)*kmult(3)+indices(6))+1 525 526 ikpt = inttoik(curicoord) 527 g0 = allg0(:,curicoord) 528 529 end subroutine get_kpt_from_indices_dense
m_double_grid/k_neighbors [ Functions ]
[ Top ] [ m_double_grid ] [ Functions ]
NAME
k_neighbors
FUNCTION
find 8 neighbors of given k-point on a coarse grid, and return them along with relative k-shift within coarse grid cell
INPUTS
kpt = k-point to be interpolated to, in full BZ kptrlatt = lattice vectors for coarse k-grid invrankkpt = rank list to find k-points
OUTPUT
rel_kpt = k-point coordinates renormalized to coarse grid cell kpt_phon_indices = indices of k-points on corners of cell
TODO
This routine is not used anymore. Deprecate or Remove?
SOURCE
888 subroutine k_neighbors (kpt, kptrlatt,krank, rel_kpt, kpt_phon_indices) 889 890 ! inputs 891 real(dp), intent(in) :: kpt(3) 892 integer, intent(in) :: kptrlatt(3,3) 893 type(krank_t), intent(in) :: krank 894 895 ! outputs 896 real(dp), intent(out) :: rel_kpt(3) 897 integer, intent(out) :: kpt_phon_indices(8) 898 899 ! local vars 900 integer :: symrankkpt 901 integer :: ir1,ir2,ir3, pr1,pr2,pr3 902 real(dp) :: redkpt(3), cornerkpt(3), res 903 904 ! ************************************************************************* 905 906 !wrap fine kpt to [0,1] 907 call wrap2_zero_one(kpt(1),redkpt(1),res) 908 call wrap2_zero_one(kpt(2),redkpt(2),res) 909 call wrap2_zero_one(kpt(3),redkpt(3),res) 910 !find 8 indices of points neighboring ikpt_phon, for interpolation 911 call interpol3d_indices (redkpt,kptrlatt(1,1),kptrlatt(2,2),kptrlatt(3,3), & 912 & ir1,ir2,ir3, pr1,pr2,pr3) 913 914 !transpose ir pr to ikpt_phon indices 915 !order of kpt_phons: 916 !ir1 ir2 ir3 917 cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/) 918 symrankkpt = krank%get_rank(cornerkpt) 919 kpt_phon_indices(1) = krank%invrank(symrankkpt) 920 !pr1 ir2 ir3 921 cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/) 922 symrankkpt = krank%get_rank (cornerkpt) 923 kpt_phon_indices(2) = krank%invrank(symrankkpt) 924 !ir1 pr2 ir3 925 cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/) 926 symrankkpt = krank%get_rank (cornerkpt) 927 kpt_phon_indices(3) = krank%invrank(symrankkpt) 928 !pr1 pr2 ir3 929 cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(ir3-1)/kptrlatt(3,3)/) 930 symrankkpt = krank%get_rank (cornerkpt) 931 kpt_phon_indices(4) = krank%invrank(symrankkpt) 932 !ir1 ir2 pr3 933 cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/) 934 symrankkpt = krank%get_rank (cornerkpt) 935 kpt_phon_indices(5) = krank%invrank(symrankkpt) 936 !pr1 ir2 pr3 937 cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(ir2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/) 938 symrankkpt = krank%get_rank (cornerkpt) 939 kpt_phon_indices(6) = krank%invrank(symrankkpt) 940 !ir1 pr2 pr3 941 cornerkpt = (/real(ir1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/) 942 symrankkpt = krank%get_rank (cornerkpt) 943 kpt_phon_indices(7) = krank%invrank(symrankkpt) 944 !pr1 pr2 pr3 945 cornerkpt = (/real(pr1-1)/kptrlatt(1,1),real(pr2-1)/kptrlatt(2,2), real(pr3-1)/kptrlatt(3,3)/) 946 symrankkpt = krank%get_rank (cornerkpt) 947 kpt_phon_indices(8) = krank%invrank(symrankkpt) 948 949 !retrieve the gkq matrix for all q, at the neighbor k vectors 950 rel_kpt(1) = redkpt(1)*kptrlatt(1,1)-real(ir1-1) 951 rel_kpt(2) = redkpt(2)*kptrlatt(2,2)-real(ir2-1) 952 rel_kpt(3) = redkpt(3)*kptrlatt(3,3)-real(ir3-1) 953 954 end subroutine k_neighbors
m_double_grid/kptfine_av [ Functions ]
[ Top ] [ m_double_grid ] [ Functions ]
NAME
kptfine_av
FUNCTION
Find the k-points of a fine grid that are around a k-point of a coarse mesh.
INPUTS
center(3) = the point of the coarse mesh around which you want know which k-points of the fine mesh belong to. qptrlatt(3,3) = qptrlatt of the considered calculation (this is obtained from the input variable ngqpt and shiftq. kpt_fine(3,nkpt_fine) = this table contain all the k-points of the fine grid in the full BZ (no sym op. allowed) and is read from the header of the dense WF file. nkpt_fine = number of k-points of the fine grid read from the header of the dense WF file.
OUTPUT
kpt_fine_sub(nkpt_sub) = k-points of the fine grid that are around center(3) nkpt_sub = number of k-points of the fine grid that are around center(3) wgt_sub(nkpt_sub) = weight of the k-points of the fine grid that are around center(3).
SOURCE
740 subroutine kptfine_av(center,qptrlatt,kpt_fine,nkpt_fine,kpt_fine_sub,nkpt_sub,wgt_sub) 741 742 !Arguments ------------------------------------ 743 !scalars 744 integer,intent(in) :: nkpt_fine 745 integer,intent(out) :: nkpt_sub 746 !arrays 747 integer,intent(in) :: qptrlatt(3,3) 748 real(dp),intent(in) :: kpt_fine(3,nkpt_fine) 749 real(dp),intent(in) :: center(3) 750 integer,pointer :: kpt_fine_sub(:) 751 real(dp),pointer :: wgt_sub(:) 752 753 !Local variables------------------------------- 754 !scalars 755 integer :: ikpt,aa,bb,cc 756 integer :: ii,jj 757 !arrays 758 real(dp) :: center_ref(3) 759 real(dp) :: kpt_fine_ref(3) 760 real(dp) :: kpt_tmp(3),kpt_tmp2(3) 761 integer,allocatable :: kpt_fine_sub_tmp(:) 762 real(dp),allocatable :: wgt_sub_tmp(:) 763 logical :: found(3) 764 765 ! ************************************************************************* 766 767 ABI_MALLOC(kpt_fine_sub_tmp,(nkpt_fine)) 768 ABI_MALLOC(wgt_sub_tmp,(nkpt_fine)) 769 770 !It is easier to work in real space using the qptrlatt matrices because in this 771 !referential any k-points sampling will be cast into an orthorhombic shape. 772 !In that space we can simply take all k-points of the fine grid that between 773 !center_ref-0.5 and center_ref+0.5 774 775 center_ref = MATMUL(qptrlatt,center) 776 777 !When considering points center(3) that lying close or on a BZ edge we need to 778 !take the k-points of the fine grid taking into account unklamp vectors. This 779 !is done with the aa, bb and cc loops. 780 781 ii = 1 782 do ikpt=1,nkpt_fine 783 kpt_tmp = kpt_fine(:,ikpt) 784 do aa=-1,1 785 kpt_tmp2(1) = kpt_tmp(1)+aa 786 do bb=-1,1 787 kpt_tmp2(2) = kpt_tmp(2)+bb 788 do cc=-1,1 789 kpt_tmp2(3) = kpt_tmp(3)+cc 790 kpt_fine_ref = MATMUL(qptrlatt,kpt_tmp2) 791 if((kpt_fine_ref(1)>=center_ref(1)-0.5-tol8).and.& 792 & (kpt_fine_ref(1)<=center_ref(1)+0.5+tol8)) then 793 if((kpt_fine_ref(2)>=center_ref(2)-0.5-tol8).and.& 794 & (kpt_fine_ref(2)<=center_ref(2)+0.5+tol8)) then 795 if((kpt_fine_ref(3)>=center_ref(3)-0.5-tol8).and.& 796 & (kpt_fine_ref(3)<=center_ref(3)+0.5+tol8)) then 797 kpt_fine_sub_tmp(ii) = ikpt 798 ii = ii +1 799 end if 800 end if 801 end if 802 end do 803 end do 804 end do 805 end do 806 807 nkpt_sub = ii-1 808 ABI_MALLOC(kpt_fine_sub,(nkpt_sub)) 809 ABI_MALLOC(wgt_sub,(nkpt_sub)) 810 811 do jj=1,nkpt_sub 812 kpt_fine_sub(jj) = kpt_fine_sub_tmp(jj) 813 end do 814 815 !We then compute a weight function. This weight function is simply a 816 !rectangular weight function that take the value 1 for k-points of the fine 817 !grid inside the cube, 0.5 for k-points that are lying on one face of the cube, 818 !0.25 for k-points that are lying on an edge of the cube and 0.125 for k-points 819 !that are lying on a peak of the cube. 820 821 wgt_sub(:) = 1.0 822 823 do ikpt=1,nkpt_sub 824 found(:) = .True. 825 kpt_tmp = kpt_fine(:,kpt_fine_sub(ikpt)) 826 do aa=-1,1 827 kpt_tmp2(1) = kpt_tmp(1)+aa 828 do bb=-1,1 829 kpt_tmp2(2) = kpt_tmp(2)+bb 830 do cc=-1,1 831 kpt_tmp2(3) = kpt_tmp(3)+cc 832 kpt_fine_ref = MATMUL(qptrlatt,kpt_tmp2) 833 if((ABS(kpt_fine_ref(1)-center_ref(1)-0.5)< tol8) .or.& 834 & (ABS(kpt_fine_ref(1)-center_ref(1)+0.5) < tol8)) then 835 if(found(1)) then 836 wgt_sub(ikpt) = wgt_sub(ikpt)*0.5 837 found(1) = .False. 838 end if 839 end if 840 if((ABS(kpt_fine_ref(2)-center_ref(2)-0.5) < tol8) .or.& 841 & (ABS(kpt_fine_ref(2)-center_ref(2)+0.5) < tol8)) then 842 if(found(2)) then 843 wgt_sub(ikpt) = wgt_sub(ikpt)*0.5 844 found(2) = .False. 845 end if 846 end if 847 if((ABS(kpt_fine_ref(3)-center_ref(3)-0.5)< tol8) .or.& 848 & (ABS(kpt_fine_ref(3)-center_ref(3)+0.5) < tol8)) then 849 if(found(3)) then 850 wgt_sub(ikpt) = wgt_sub(ikpt)*0.5 851 found(3) = .False. 852 end if 853 end if 854 end do 855 end do 856 end do 857 end do 858 859 ABI_FREE(kpt_fine_sub_tmp) 860 ABI_FREE(wgt_sub_tmp) 861 862 end subroutine kptfine_av