1! Copyright (C) 2017 Quantum ESPRESSO Foundation 2! Author: Ivan Carnimeo 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!-------------------------------------- 9MODULE loc_scdm_k 10 !-------------------------------------- 11 !! Variables and subroutines for localizing molecular orbitals based 12 !! on a modified SCDM approach. Original SCDM method: 13 ! 14 !! A. Damle, L. Lin, L. Ying: J. Comput. Phys., 334 (2017) 1-5 15 ! 16 !! Ivan Carnimeo: SCDM has been implemented with a parallel prescreening algorithm 17 !! in order to save CPU time and memory for large grids. 18 ! 19 USE kinds, ONLY : DP 20 USE io_global, ONLY : stdout 21 USE exx, ONLY : dfftt, exxbuff, exxmat 22 USE exx_base, ONLY : nkqs 23 ! 24 IMPLICIT NONE 25 ! 26 SAVE 27 ! 28 REAL(DP), PARAMETER :: Zero=0.0_DP, One=1.0_DP, Two=2.0_DP, Three=2.0_DP 29 INTEGER :: n_scdm = 1 30 ! 31 CONTAINS 32 ! 33!------------------------------------------------------------------------ 34SUBROUTINE localize_orbitals_k( ) 35 !---------------------------------------------------------------------- 36 !! Driver for SCDM_k orbital localization: 37 ! 38 !! * \(\texttt{SCDM_PGG_k}\): performs SCDM localization using the 39 !! parallel prescreening algorithm; 40 !! * \(\texttt{measure_localization_k}\): computes the absolute overlap 41 !! integrals required by \(\texttt{vexx_loc_k in}\) 42 !! exx.f90 and puts them in \(\texttt{exxmat}\). 43 ! 44 !! NOTE: \(\texttt{localize_orbitals_k}\) does not have iscdm/nscdm because 45 !! alignment is not implemented for K-Points. 46 ! 47 USE noncollin_module, ONLY : npol 48 USE exx, ONLY : x_occupation 49 USE klist, ONLY : nks 50 USE wvfct, ONLY : nbnd 51 USE mp_bands, ONLY : me_bgrp 52 USE control_flags, ONLY : lmd 53 ! 54 IMPLICIT NONE 55 ! 56 INTEGER :: NGrid, ikq, jk, k, NBands, nnk 57 CHARACTER(LEN=1) :: HowTo 58 REAL(DP) :: Tot1, Tot2, Ave1, tot, ave 59 ! 60 IF ( n_scdm /= 1 ) CALL errore( 'localize_orbitals_k','nscdm for K-points NYI.', 1 ) 61 IF ( lmd ) CALL errore( 'localize_orbitals_k', 'localization with K-points not tested.', 1 ) 62 ! 63 NGrid = dfftt%nnr * npol 64 HowTo = 'G' ! How to compute the absolute overlap integrals (only G for k) 65 ! 66 exxmat = One ! initialized to one because the absolute overlaps between 67 ! occupied and virtual orbitals are assumed always large 68 ! and all ov integrals will be computed by vexx_loc_k. 69 ! 70 NBands = INT(SUM(x_occupation(:,1))) 71 ! 72 WRITE (stdout,*) ' ' 73 WRITE (stdout,*) 'NBands = ', NBands, ' nks = ', nks, ' nkqs = ', nkqs 74 WRITE (stdout,'(5X,A)') 'Canonical Orbitals ' 75 ! 76 Tot1 = Zero 77 Ave1 = Zero 78 Tot2 = Zero 79 nnk = 0 80 ! 81 DO ikq = 1, nkqs 82 ! 83 CALL measure_localization_k( NBands, ikq, tot, ave ) 84 ! 85 Tot1 = Tot1 + tot 86 Ave1 = Ave1 + ave 87 ! 88 DO jk = 1, nks 89 ! 90 nnk = nnk + 1 91 ! 92 CALL AbsOvG_k( NBands, ikq, jk, tot, ave ) 93 ! 94 Tot2 = Tot2 + tot 95 ! 96 ENDDO 97 ! 98 ENDDO 99 ! 100 Ave1 = Ave1/DBLE(nkqs) 101 WRITE (stdout,'(7X,A,f24.6)') 'Total AbsOv =', Tot2 102 WRITE (stdout,'(7X,A,f24.6)') 'Aver. AbsOv =', Tot2/DBLE(nnk) 103 WRITE (stdout,'(7X,A,f24.6)') 'Total Spread [A**2] =', Tot1 104 WRITE (stdout,'(7X,A,f24.6)') 'Aver. Spread [A**2] =', Ave1 105 ! 106 ! IF (me_bgrp == 0) THEN 107 ! CALL AbsOv_histogram_k( NBands, 'hist1.dat' ) 108 ! ENDIF 109 ! 110 WRITE (stdout,'(5X,A)') 'SCDM-PGG_k localization' 111 DO ikq = 1, nkqs 112 CALL SCDM_PGG_k( NGrid, NBands, ikq ) 113 ENDDO 114 WRITE (stdout,'(7X,A)') 'SCDM-PGG_k done ' 115 ! 116 WRITE (stdout,'(5X,A)') 'Localized Orbitals ' 117 ! 118 Tot1 = Zero 119 Ave1 = Zero 120 Tot2 = Zero 121 nnk = 0 122 ! 123 DO ikq = 1, nkqs 124 ! 125 CALL measure_localization_k( NBands, ikq, tot, ave ) 126 ! 127 Tot1 = Tot1 + tot 128 Ave1 = Ave1 + ave 129 ! 130 DO jk = 1, nks 131 ! 132 nnk = nnk + 1 133 ! 134 CALL AbsOvG_k( NBands, ikq, jk, tot, ave ) 135 ! 136 Tot2 = Tot2 + tot 137 ! 138 ENDDO 139 ! 140 ENDDO 141 ! 142 Ave1 = Ave1 / DBLE(nkqs) 143 WRITE (stdout,'(7X,A,f24.6)') 'Total AbsOv =', Tot2 144 WRITE (stdout,'(7X,A,f24.6)') 'Aver. AbsOv =', Tot2/DBLE(nnk) 145 WRITE (stdout,'(7X,A,f24.6)') 'Total Spread [A**2] =', Tot1 146 WRITE (stdout,'(7X,A,f24.6)') 'Aver. Spread [A**2] =', Ave1 147 ! 148 ! IF (me_bgrp == 0) THEN 149 ! CALL AbsOv_histogram_k( NBands, 'hist2.dat' ) 150 ! ENDIF 151 ! 152END SUBROUTINE localize_orbitals_k 153! 154! 155!------------------------------------------------------------------------------- 156SUBROUTINE measure_localization_k( NBands, IK, TotSpread, AveSpread ) 157 !----------------------------------------------------------------------------- 158 !! Computes the absolute overlap integrals required by \(\texttt{vexx_loc_k}\) 159 !! in exx.f90 and puts them in \(\texttt{exxmat}\). 160 ! 161 USE noncollin_module, ONLY : npol 162 USE cell_base, ONLY : alat, omega, at, bg 163 USE exx, ONLY : compute_density_k 164 USE constants, ONLY : bohr_radius_angs 165 ! 166 IMPLICIT NONE 167 ! 168 REAL(DP), INTENT(OUT) :: TotSpread 169 !! Total Spread increment per IK 170 REAL(DP), INTENT(OUT) :: AveSpread 171 !! Average Spread increment per IK 172 INTEGER :: NBands 173 !! number of bands (with auxiliary functions) 174 INTEGER :: IK 175 !! index on k+q points 176 ! 177 ! ... local variables 178 ! 179 INTEGER :: ibnd 180 REAL(DP) :: RJunk(4), SpreadPBC(3) 181 ! 182 CALL start_clock( 'measure' ) 183 ! 184 TotSpread = Zero 185 AveSpread = Zero 186 ! 187 DO ibnd = 1, NBands 188 ! 189 CALL compute_density_k( .FALSE., .FALSE., RJunk(1:3), SpreadPBC, RJunk(4), & 190 exxbuff(1,ibnd,IK), exxbuff( 1,ibnd,IK), & 191 dfftt%nnr*npol, ibnd, ibnd ) 192 TotSpread = TotSpread + SpreadPBC(1) + SpreadPBC(2) + SpreadPBC(3) 193 ! 194 ENDDO 195 ! 196 TotSpread = TotSpread * bohr_radius_angs**2 197 AveSpread = TotSpread / DBLE(NBands) 198 ! 199 ! WRITE (stdout,'(7X,A,I3)') 'IK = ', IK 200 ! WRITE (stdout,'(7X,A,f12.6,I3)') 'Total Spread [A**2] =', TotSpread 201 ! WRITE (stdout,'(7X,A,f12.6,I3)') 'Aver. Spread [A**2] =', AveSpread 202 ! 203 CALL stop_clock( 'measure' ) 204 ! 205END SUBROUTINE measure_localization_k 206! 207! 208!--------------------------------------------------------------------- 209SUBROUTINE AbsOvG_k( NBands, IKQ, JK, loc_diag, loc_off ) 210 !--------------------------`---------------------------------------- 211 !! Computes the Absolute Overlap in G-space (cutoff might not be 212 !! accurate for the moduli of the wavefunctions). 213 ! 214 USE noncollin_module, ONLY : npol 215 USE fft_interfaces, ONLY : fwfft 216 USE wvfct, ONLY : npwx 217 USE exx_band, ONLY : igk_exx 218 USE klist, ONLY : nks, ngk 219 ! 220 IMPLICIT NONE 221 ! 222 INTEGER :: NBands 223 !! number of bands (with auxiliary functions) 224 INTEGER :: IKQ 225 !! index on k+q points 226 INTEGER :: JK 227 !! index on k-point 228 REAL(DP), INTENT(OUT) :: loc_diag 229 !! Total Charge 230 REAL(DP), INTENT(OUT) :: loc_off 231 !! Total Abs. Overlap 232 ! 233 ! ... local variables 234 ! 235 INTEGER :: ibnd, jbnd, ig, kk, npw 236 REAL(DP) :: tmp 237 COMPLEX(DP), ALLOCATABLE :: buffer(:), GorbtI(:,:), GorbtJ(:,:) 238 INTEGER, EXTERNAL :: global_kpoint_index 239 COMPLEX(DP), ALLOCATABLE :: Mat(:,:) 240 ! 241 CALL start_clock( 'measure' ) 242 ! 243 ! WRITE (stdout, '(5X,A)') ' ' 244 ! WRITE (stdout, '(5X,A)') 'Absolute Overlap calculated in G-space (K)' 245 ! 246 kk = global_kpoint_index ( nks, JK ) 247 ! 248 ALLOCATE( Mat(NBands,NBands) ) 249 ALLOCATE( buffer(dfftt%nnr*npol), GorbtI(npwx,NBands), GorbtJ(npwx,NBands) ) 250 ! 251 GorbtJ = (Zero,Zero) 252 GorbtI = (Zero,Zero) 253 npw = ngk(JK) 254 ! 255 DO jbnd = 1, NBands 256 ! 257 buffer(:) = ABS(exxbuff(:,jbnd,JK)) 258 ! buffer(:) = exxbuff(:,jbnd,JK) 259 ! 260 CALL fwfft( 'Wave' , buffer, dfftt ) 261 ! 262 DO ig = 1, npw 263 GorbtJ(ig,jbnd) = buffer(dfftt%nl(igk_exx(ig,kk))) 264 ENDDO 265 ! 266 buffer(:) = ABS(exxbuff(:,jbnd,IKQ)) 267 ! buffer(:) = exxbuff(:,jbnd,IKQ) 268 ! 269 CALL fwfft( 'Wave' , buffer, dfftt ) 270 ! 271 DO ig = 1, npw 272 GorbtI(ig,jbnd) = buffer(dfftt%nl(igk_exx(ig,kk))) 273 ENDDO 274 ! 275 ENDDO 276 ! 277 ! IF (IKQ == JK) THEN 278 ! CALL matcalc_k( 'AbsOv-',.FALSE., 2, 0, npwx, NBands, NBands, GorbtI, GorbtJ, Mat, tmp ) 279 ! ELSE 280 CALL matcalc_k( 'AbsOv-',.FALSE., 0, 0, npwx, NBands, NBands, GorbtI, GorbtJ, Mat, tmp ) 281 ! END IF 282 ! 283 loc_diag = Zero 284 loc_off = Zero 285 ! 286 DO ibnd = 1, NBands 287 ! 288 loc_diag = loc_diag + DBLE(Mat(ibnd,ibnd)) 289 ! 290 DO jbnd = 1, ibnd-1 291 ! 292 loc_off = loc_off + DBLE(Mat(ibnd,jbnd)) + DBLE(Mat(jbnd,ibnd)) 293 exxmat(ibnd, IKQ, jbnd, JK) = DBLE(Mat(ibnd,jbnd)) 294 exxmat(jbnd, IKQ, ibnd, JK) = DBLE(Mat(jbnd,ibnd)) 295 ! 296 ENDDO 297 ! 298 ENDDO 299 ! 300 WRITE (stdout,'(7X,5(A,I3),2(A,f12.6))') 'IKQ = ', IKQ, & 301 ' JK = ', JK, ' kk = ', kk, ' NBands = ', NBands, ' size = ', SIZE(exxbuff,3), & 302 ' Total Charge =', loc_diag, ' Total Abs. Overlap =', loc_off 303 ! 304 DEALLOCATE( buffer, GorbtI, GorbtJ ) 305 DEALLOCATE( Mat ) 306 ! 307 CALL stop_clock( 'measure' ) 308 ! 309END SUBROUTINE AbsOvG_k 310! 311! 312!------------------------------------------------------------------------------------ 313SUBROUTINE AbsOv_histogram_k( n, filename ) 314 !---------------------------------------------------------------------------------- 315 !! Prints the distribution of absolute overlap integrals (as an histogram) in the 316 !! range [0,1] with NHist intervals. 317 ! 318 !! WARNING: HUGE amount of integrals, use only for VERY small systems with a few 319 !! K-Points. 320 ! 321 USE klist, ONLY: nks 322 ! 323 IMPLICIT NONE 324 ! 325 CHARACTER(LEN=*), INTENT(IN) :: filename 326 !! file with histogram data 327 INTEGER, INTENT(IN) :: n 328 !! number of bands 329 ! 330 ! ... local variables 331 ! 332 INTEGER :: i,j,k, iik, jjk, NHist 333 INTEGER, ALLOCATABLE :: histogram(:) 334 REAL(DP), ALLOCATABLE :: XHist(:) 335 REAL(DP) :: xstart, xstep, integral 336 INTEGER :: io_histogram 337 INTEGER, EXTERNAL :: find_free_unit 338 ! 339 NHist = 1000 340 xstep = One/FLOAT(NHist) 341 xstart = One/FLOAT(NHist)/Two 342 ! 343 WRITE (stdout,'(A,I7,2(A,f12.6))') 'NHist = ', NHist , ' xstep = ', xstep, ' xstart = ', xstart 344 ! 345 ALLOCATE( histogram(NHist), XHist(NHist) ) 346 ! 347 XHist = Zero 348 histogram = 0 349 ! 350 DO k = 1, NHist 351 XHist(k) = xstart + FLOAT(k-1)*xstep 352 ENDDO 353 ! 354 DO i = 1, n 355 DO iik = 1, nkqs 356 ! 357 DO j = 1, n 358 DO jjk = 1, nks 359 ! 360 integral = exxmat(i, iik, j, jjk) 361 ! 362 IF (integral < Zero) THEN 363 CALL errore( 'AbsOv_histogram_k','Abs. Ov. < 0 found.', 1 ) 364 ELSE 365 DO k = 1, NHist 366 IF (integral >= (XHist(k)-xstart) .AND. integral < (XHist(k)+xstart)) & 367 histogram(k) = histogram(k) + 1 368 ENDDO 369 ENDIF 370 ! 371 ENDDO 372 ENDDO 373 ! 374 ENDDO 375 ENDDO 376 ! 377 io_histogram = find_free_unit() 378 OPEN( io_histogram, file=filename, status='unknown' ) 379 ! 380 DO k = 1, NHist 381 WRITE (io_histogram, '(f12.6,2I10)') XHist(k), histogram(k) 382 ENDDO 383 ! 384 CLOSE(io_histogram) 385 ! 386 DEALLOCATE( histogram, XHist ) 387 ! 388END SUBROUTINE AbsOv_histogram_k 389! 390! 391!----------------------------------------------------------------------------- 392SUBROUTINE SCDM_PGG_k( NGrid, NBands, IKQ ) 393 !--------------------------------------------------------------------------- 394 !! Density matrix localization (I/O in psi). K-Points version. 395 ! 396 USE mp_bands, ONLY: nproc_bgrp 397 USE loc_scdm, ONLY: scdm_thresholds, scdm_points, scdm_prescreening 398 ! 399 IMPLICIT NONE 400 ! 401 INTEGER, INTENT(IN) :: NGrid 402 !! Number of grid points 403 INTEGER, INTENT(IN) :: NBands 404 !! Number of bands 405 INTEGER, INTENT(IN) :: IKQ 406 !! Index on k+q points 407 ! 408 ! ... local variables 409 ! 410 COMPLEX(DP), ALLOCATABLE :: QRbuff(:,:), mat(:,:) 411 INTEGER, ALLOCATABLE :: pivot(:), list(:), cpu_npt(:) 412 REAL(DP), ALLOCATABLE :: den(:), grad_den(:,:) 413 INTEGER :: nptot 414 REAL(DP) :: ThrDen, ThrGrd 415 ! 416 CALL start_clock( 'localization' ) 417 ! 418 ! WRITE (stdout,'(5X,A)') 'SCDM localization with prescreening' 419 ! 420 ALLOCATE( den(dfftt%nnr), grad_den(3, dfftt%nnr) ) 421 ! 422 CALL scdm_thresholds( den, grad_den, ThrDen, ThrGrd ) 423 ! 424 ALLOCATE( cpu_npt(0:nproc_bgrp-1) ) 425 ! 426 CALL scdm_points( den, grad_den, ThrDen, ThrGrd, cpu_npt, nptot ) 427 ! 428 ALLOCATE( list(nptot), pivot(nptot) ) 429 ! 430 CALL scdm_prescreening_k( NGrid, NBands, exxbuff(1,1,IKQ), den, grad_den, & 431 ThrDen, ThrGrd, cpu_npt, nptot, list, pivot ) 432 DEALLOCATE( den, grad_den ) 433 ! 434 ! Psi(pivot(1:NBands),:) in mat 435 ! 436 ALLOCATE( mat(NBands,NBands) ) 437 ! 438 mat = (Zero, Zero) 439 ! 440 CALL scdm_fill_k( .TRUE., nptot, NGrid, NBands, cpu_npt, pivot, list, & 441 exxbuff(1,1,IKQ), Mat ) 442 ! 443 ! Pc = Psi * Psi(pivot(1:NBands),:)' in QRbuff 444 ! 445 ALLOCATE( QRbuff(NGrid, NBands) ) 446 ! 447 QRbuff = (Zero, Zero) 448 ! 449 CALL ZGEMM( 'N' , 'N' , NGrid, NBands, NBands, (One,Zero), exxbuff(1,1,IKQ), & 450 NGrid, mat, NBands, (Zero,Zero), QRbuff, NGrid ) 451 ! 452 ! Orthonormalization 453 ! Pc(pivot(1:NBands),:) in mat 454 ! 455 mat = (Zero, Zero) 456 ! 457 CALL scdm_fill_k( .FALSE., nptot, NGrid, NBands, cpu_npt, pivot, list, QRBuff, mat ) 458 ! 459 DEALLOCATE( cpu_npt ) 460 ! 461 ! Cholesky(psi)^(-1) in mat 462 CALL invchol_k(NBands,mat) 463 ! 464 ! Phi = Pc * Chol^(-1) = QRbuff * mat 465 CALL ZGEMM( 'N' , 'T' , NGrid, NBands, NBands, (One,Zero), QRbuff, NGrid, mat, & 466 NBands, (Zero,Zero), exxbuff(1,1,IKQ), NGrid ) 467 ! 468 DEALLOCATE( QRbuff, mat, pivot, list ) 469 ! 470 ! WRITE (stdout,'(7X,A)') 'SCDM-PGG_k done ' 471 ! 472 CALL stop_clock( 'localization' ) 473 ! 474END SUBROUTINE SCDM_PGG_k 475! 476! 477!----------------------------------------------------------------------------------------------------- 478SUBROUTINE scdm_prescreening_k( NGrid, NBands, psi, den, grad_den, ThrDen, ThrGrd, & 479 cpu_npt, nptot, list, pivot ) 480 !--------------------------------------------------------------------------------------------------- 481 !! Get List from ThrDen and ThrGrd, and Pivot from the QRCP of small. 482 ! 483 USE mp, ONLY : mp_sum 484 USE mp_bands, ONLY : intra_bgrp_comm, me_bgrp, nproc_bgrp 485 ! 486 IMPLICIT NONE 487 ! 488 INTEGER, INTENT(OUT) :: list(nptot) 489 INTEGER, INTENT(OUT) :: pivot(nptot) 490 INTEGER, INTENT(IN) :: cpu_npt(0:nproc_bgrp-1) 491 !! Number of relevant points per processor 492 INTEGER, INTENT(IN) :: nptot 493 !! Number of relevant points for the allocation 494 INTEGER, INTENT(IN) :: NGrid 495 !! Number of grid points 496 INTEGER, INTENT(IN) :: NBands 497 !! Number of bands 498 COMPLEX(DP), INTENT(IN) :: psi(NGrid,NBands) 499 REAL(DP), INTENT(IN) :: den(dfftt%nnr) 500 REAL(DP), INTENT(IN) :: grad_den(3,dfftt%nnr) 501 REAL(DP), INTENT(IN) :: ThrDen 502 REAL(DP), INTENT(IN) :: ThrGrd 503 ! 504 ! ... local variables 505 ! 506 INTEGER :: ir, ir_end, INFO, lwork 507 INTEGER :: n, npt, ncpu_start 508 REAL(DP) :: grad 509 COMPLEX(DP), ALLOCATABLE :: small(:,:), tau(:), work(:) 510 REAL(DP), ALLOCATABLE :: rwork(:) 511 ! 512#if defined (__MPI) 513 ir_end = dfftt%nr1x*dfftt%my_nr2p*dfftt%my_nr3p 514#else 515 ir_end = dfftt%nnr 516#endif 517 ! 518 ! find the map of the indeces 519 ALLOCATE( small(NBands,nptot) ) 520 small = (Zero, Zero) 521 list = 0 522 n = 0 523 ! 524 DO ir = 1, ir_end 525 grad = One 526 IF (den(ir) > ThrDen) THEN 527 grad = SQRT( grad_den(1,ir)**2 + grad_den(2,ir)**2 + grad_den(3,ir)**2 ) 528 IF (grad < ThrGrd) THEN 529 n = n + 1 530 ncpu_start = SUM(cpu_npt(0:me_bgrp-1)) 531 npt = ncpu_start+n 532 small(:,npt) = psi(ir,:) 533 list(npt) = ir 534 ENDIF 535 ENDIF 536 ENDDO 537 ! 538 CALL mp_sum( small, intra_bgrp_comm ) 539 CALL mp_sum( list, intra_bgrp_comm ) 540 ! 541 ! perform the QRCP on the small matrix and get pivot 542 lwork = 4*nptot 543 ALLOCATE( tau(nptot), work(lwork), rwork(2*nptot) ) 544 tau = (Zero, Zero) 545 work = (Zero, Zero) 546 pivot = 0 547 INFO = -1 548 CALL ZGEQP3( NBands, nptot, small, NBands, pivot, tau, work, lwork, rwork, INFO ) 549 DEALLOCATE( tau, work, rwork, small ) 550 ! 551END SUBROUTINE scdm_prescreening_k 552! 553! 554!---------------------------------------------------------------------------------- 555SUBROUTINE scdm_fill_k( op, nptot, NGrid, NBands, CPUPts, Pivot, List, Vect, Mat ) 556 !-------------------------------------------------------------------------------- 557 !! Fills the matrix Mat with the elements of Vect mapped by CPUPts, Pivot 558 !! and List. 559 ! 560 USE mp, ONLY: mp_sum 561 USE mp_bands, ONLY: intra_bgrp_comm, me_bgrp, nproc_bgrp 562 ! 563 IMPLICIT NONE 564 ! 565 INTEGER, INTENT(IN) :: NBands 566 INTEGER, INTENT(IN) :: NGrid 567 INTEGER, INTENT(IN) :: nptot 568 INTEGER, INTENT(IN) :: CPUPts(0:nproc_bgrp-1) 569 INTEGER, INTENT(IN) :: Pivot(nptot) 570 INTEGER, INTENT(IN) :: List(nptot) 571 COMPLEX(DP), INTENT(IN) :: Vect(NGrid,NBands) 572 COMPLEX(DP), INTENT(OUT) :: Mat(NBands,NBands) 573 ! 574 ! ... local variables 575 ! 576 INTEGER :: i, NStart, NEnd 577 LOGICAL :: op 578 ! 579 Mat = (Zero, Zero) 580 ! 581 DO i = 1, NBands 582 NStart = SUM(CPUPts(0:me_bgrp-1)) 583 NEnd = SUM(CPUPts(0:me_bgrp)) 584 IF (Pivot(i) <= NEnd .AND. Pivot(i) >= NStart+1) THEN 585 IF (op) THEN 586 Mat(:,i) = CONJG(Vect(List(pivot(i)),:)) 587 ELSE 588 Mat(:,i) = Vect(List(pivot(i)),:) 589 ENDIF 590 ENDIF 591 ENDDO 592 ! 593 CALL mp_sum( Mat, intra_bgrp_comm ) 594 ! 595END SUBROUTINE scdm_fill_k 596! 597! 598! 599END MODULE loc_scdm_k 600!----------------------------------------------------------------------- 601