1!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 2! Copyright 2010. Los Alamos National Security, LLC. This material was ! 3! produced under U.S. Government contract DE-AC52-06NA25396 for Los Alamos ! 4! National Laboratory (LANL), which is operated by Los Alamos National ! 5! Security, LLC for the U.S. Department of Energy. The U.S. Government has ! 6! rights to use, reproduce, and distribute this software. NEITHER THE ! 7! GOVERNMENT NOR LOS ALAMOS NATIONAL SECURITY, LLC MAKES ANY WARRANTY, ! 8! EXPRESS OR IMPLIED, OR ASSUMES ANY LIABILITY FOR THE USE OF THIS ! 9! SOFTWARE. If software is modified to produce derivative works, such ! 10! modified software should be clearly marked, so as not to confuse it ! 11! with the version available from LANL. ! 12! ! 13! Additionally, this program is free software; you can redistribute it ! 14! and/or modify it under the terms of the GNU General Public License as ! 15! published by the Free Software Foundation; version 2.0 of the License. ! 16! Accordingly, this program is distributed in the hope that it will be ! 17! useful, but WITHOUT ANY WARRANTY; without even the implied warranty of ! 18! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General ! 19! Public License for more details. ! 20!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 21 22SUBROUTINE KBLDNEWH 23 24 USE CONSTANTS_MOD 25 USE SETUPARRAY 26 USE NEBLISTARRAY 27 USE XBOARRAY 28 USE NONOARRAY 29 USE UNIVARRAY 30 USE KSPACEARRAY 31 USE MYPRECISION 32 33 IMPLICIT NONE 34 35 INTEGER :: I, J, NEWJ, K, L, II, JJ, KK, MM, MP, NN, SUBI 36 INTEGER :: IBRA, IKET, LBRA, LKET, MBRA, MKET 37 INTEGER :: INDEX, INDI, INDJ 38 INTEGER :: SWITCH, PREVJ 39 INTEGER :: PBCI, PBCJ, PBCK 40 INTEGER :: BASISI(5), BASISJ(5), LBRAINC, LKETINC 41 INTEGER :: NORBI, KCOUNT 42 INTEGER :: KX, KY, KZ 43 REAL(LATTEPREC) :: KX0, KY0, KZ0 44 REAL(LATTEPREC) :: ALPHA, BETA, COSBETA, PHI, TMP, PERM 45 REAL(LATTEPREC) :: RIJ(3), MAGR2, MAGR, MAGRP, RCUTTB 46 REAL(LATTEPREC) :: MAXARRAY(20), MAXRCUT, MAXRCUT2 47 REAL(LATTEPREC) :: ANGFACTOR 48 REAL(LATTEPREC) :: KPOINT(3), KDOTL 49 REAL(LATTEPREC) :: AMMBRA, WIGLBRAMBRA 50 COMPLEX(LATTEPREC) :: BLOCH, KHTMP 51 REAL(LATTEPREC), EXTERNAL :: UNIVSCALE, WIGNERD, SLMMP, TLMMP, AM, BM 52 53 HK = CMPLX(ZERO) 54 55 INDEX = 0 56 57 ! Build diagonal elements 58 DO I = 1, NATS 59 60 SELECT CASE(BASIS(ELEMPOINTER(I))) 61 62 CASE("s") 63 64 INDEX = INDEX + 1 65 HK(INDEX, INDEX, 1) = CMPLX(HES(ELEMPOINTER(I))) 66 67 CASE("p") 68 69 DO SUBI = 1, 3 70 INDEX = INDEX + 1 71 HK(INDEX,INDEX, 1) = CMPLX(HEP(ELEMPOINTER(I))) 72 ENDDO 73 74 CASE("d") 75 76 DO SUBI = 1, 5 77 INDEX = INDEX + 1 78 HK(INDEX,INDEX, 1) = CMPLX(HED(ELEMPOINTER(I))) 79 ENDDO 80 81 CASE("f") 82 83 DO SUBI = 1, 7 84 INDEX = INDEX + 1 85 HK(INDEX,INDEX, 1) = CMPLX(HEF(ELEMPOINTER(I))) 86 ENDDO 87 88 CASE("sp") 89 90 DO SUBI = 1, 4 91 92 INDEX = INDEX + 1 93 IF (SUBI .EQ. 1) THEN 94 HK(INDEX,INDEX, 1) = CMPLX(HES(ELEMPOINTER(I))) 95 ELSE 96 HK(INDEX,INDEX, 1) = CMPLX(HEP(ELEMPOINTER(I))) 97 ENDIF 98 99 ENDDO 100 101 CASE("sd") 102 103 DO SUBI = 1, 6 104 105 INDEX = INDEX + 1 106 IF (SUBI .EQ. 1) THEN 107 HK(INDEX,INDEX, 1) = CMPLX(HES(ELEMPOINTER(I))) 108 ELSE 109 HK(INDEX,INDEX, 1) = CMPLX(HED(ELEMPOINTER(I))) 110 ENDIF 111 112 ENDDO 113 114 CASE("sf") 115 116 DO SUBI = 1, 8 117 118 INDEX = INDEX + 1 119 IF (SUBI .EQ. 1) THEN 120 HK(INDEX,INDEX, 1) = CMPLX(HES(ELEMPOINTER(I))) 121 ELSE 122 HK(INDEX,INDEX, 1) = CMPLX(HEF(ELEMPOINTER(I))) 123 ENDIF 124 125 ENDDO 126 127 CASE("pd") 128 129 DO SUBI = 1, 8 130 131 INDEX = INDEX + 1 132 IF (SUBI .LE. 3) THEN 133 HK(INDEX,INDEX, 1) = CMPLX(HEP(ELEMPOINTER(I))) 134 ELSE 135 HK(INDEX,INDEX, 1) = CMPLX(HED(ELEMPOINTER(I))) 136 ENDIF 137 138 ENDDO 139 140 CASE("pf") 141 142 DO SUBI = 1, 10 143 144 INDEX = INDEX + 1 145 IF (SUBI .LE. 3) THEN 146 HK(INDEX,INDEX, 1) = CMPLX(HEP(ELEMPOINTER(I))) 147 ELSE 148 HK(INDEX,INDEX, 1) = CMPLX(HEF(ELEMPOINTER(I))) 149 ENDIF 150 151 ENDDO 152 153 CASE("df") 154 155 DO SUBI = 1, 12 156 157 INDEX = INDEX + 1 158 IF (SUBI .LE. 5) THEN 159 HK(INDEX,INDEX, 1) = CMPLX(HED(ELEMPOINTER(I))) 160 ELSE 161 HK(INDEX,INDEX, 1) = CMPLX(HEF(ELEMPOINTER(I))) 162 ENDIF 163 164 ENDDO 165 166 CASE("spd") 167 168 DO SUBI = 1, 9 169 170 INDEX = INDEX + 1 171 IF (SUBI .EQ. 1) THEN 172 HK(INDEX,INDEX, 1) = CMPLX(HES(ELEMPOINTER(I))) 173 ELSEIF (SUBI .GT. 1 .AND. SUBI .LE. 4) THEN 174 HK(INDEX,INDEX, 1) = CMPLX(HEP(ELEMPOINTER(I))) 175 ELSE 176 HK(INDEX,INDEX, 1) = CMPLX(HED(ELEMPOINTER(I))) 177 ENDIF 178 179 ENDDO 180 181 CASE("spf") 182 183 DO SUBI = 1, 11 184 185 INDEX = INDEX + 1 186 IF (SUBI .EQ. 1) THEN 187 HK(INDEX,INDEX, 1) = CMPLX(HES(ELEMPOINTER(I))) 188 ELSEIF (SUBI .GT. 1 .AND. SUBI .LE. 4) THEN 189 HK(INDEX,INDEX, 1) = CMPLX(HEP(ELEMPOINTER(I))) 190 ELSE 191 HK(INDEX,INDEX, 1) = CMPLX(HEF(ELEMPOINTER(I))) 192 ENDIF 193 194 ENDDO 195 196 CASE("sdf") 197 198 DO SUBI = 1, 13 199 200 INDEX = INDEX + 1 201 IF (SUBI .EQ. 1) THEN 202 HK(INDEX,INDEX, 1) = CMPLX(HES(ELEMPOINTER(I))) 203 ELSEIF (SUBI .GT. 1 .AND. SUBI .LE. 6) THEN 204 HK(INDEX,INDEX, 1) = CMPLX(HED(ELEMPOINTER(I))) 205 ELSE 206 HK(INDEX,INDEX, 1) = CMPLX(HEF(ELEMPOINTER(I))) 207 ENDIF 208 209 ENDDO 210 211 212 CASE("pdf") 213 214 DO SUBI = 1, 15 215 216 INDEX = INDEX + 1 217 IF (SUBI .LE. 3) THEN 218 HK(INDEX,INDEX, 1) = CMPLX(HEP(ELEMPOINTER(I))) 219 ELSEIF (SUBI .GT. 3 .AND. SUBI .LE. 8) THEN 220 HK(INDEX,INDEX, 1) = CMPLX(HED(ELEMPOINTER(I))) 221 ELSE 222 HK(INDEX,INDEX, 1) = CMPLX(HEF(ELEMPOINTER(I))) 223 ENDIF 224 225 ENDDO 226 227 CASE("spdf") 228 229 DO SUBI = 1, 16 230 231 INDEX = INDEX + 1 232 IF (SUBI .EQ. 1) THEN 233 HK(INDEX, INDEX, 1) = CMPLX(HES(ELEMPOINTER(I))) 234 ELSEIF (SUBI .GT. 1 .AND. SUBI .LE. 4) THEN 235 HK(INDEX, INDEX, 1) = CMPLX(HEP(ELEMPOINTER(I))) 236 ELSEIF (SUBI .GT. 1 .AND. SUBI .LE. 4) THEN 237 HK(INDEX, INDEX, 1) = CMPLX(HEP(ELEMPOINTER(I))) 238 ELSEIF (SUBI .GT. 4 .AND. SUBI .LE. 9) THEN 239 HK(INDEX, INDEX, 1) = CMPLX(HED(ELEMPOINTER(I))) 240 ELSE 241 HK(INDEX, INDEX, 1) = CMPLX(HEF(ELEMPOINTER(I))) 242 ENDIF 243 244 ENDDO 245 246 END SELECT 247 248 ENDDO 249 250 DO I = 2, NKTOT 251 DO J = 1, HDIM 252 HK(J,J,I) = HK(J,J,1) 253 ENDDO 254 ENDDO 255 256 ! We assign the diagonal elements in ADDQDEP 257 258 ! IF (BASISTYPE .EQ. "NONORTHO") THEN 259 260 ! SMAT = ZERO 261 262 ! DO I = 1, HDIM 263 ! SMAT(I,I) = ONE 264 ! ENDDO 265 266 ! ENDIF 267 268 ! Loop over Kpoints 269 270 271 KX0 = PI*(ONE - REAL(NKX))/(REAL(NKX)*BOXDIMS(1)) - PI*KSHIFT(1) 272 KY0 = PI*(ONE - REAL(NKY))/(REAL(NKY)*BOXDIMS(2)) - PI*KSHIFT(2) 273 KZ0 = PI*(ONE - REAL(NKZ))/(REAL(NKZ)*BOXDIMS(3)) - PI*KSHIFT(3) 274 275 ! KCOUNT = 0 276 277 ! DO KX = 1, NKX 278 279 ! KPOINT(1) = KX0 + TWO*PI*REAL(KX-1)/(BOXDIMS(1)*REAL(NKX)) 280 281 ! DO KY = 1, NKY 282 283 ! KPOINT(2) = KY0 + TWO*PI*REAL(KY-1)/(BOXDIMS(2)*REAL(NKY)) 284 285 ! DO KZ = 1, NKZ 286 287 ! KPOINT(3) = KZ0 + TWO*PI*REAL(KZ-1)/(BOXDIMS(3)*REAL(NKZ)) 288 289 ! KCOUNT = KCOUNT+1 290 291 !$OMP PARALLEL DO DEFAULT (NONE) & 292 !$OMP SHARED(KPOINT, KCOUNT, NATS, BASIS, ELEMPOINTER, TOTNEBTB, NEBTB) & 293 !$OMP SHARED(CR, BOXDIMS, HK, NOINT, ATELE, ELE1, ELE2) & 294 !$OMP SHARED(BOND, OVERL, MATINDLIST, BASISTYPE) & 295 !$OMP PRIVATE(I, J, K, NEWJ, BASISI, BASISJ, INDI, INDJ, PBCI, PBCJ, PBCK) & 296 !$OMP PRIVATE(RIJ, MAGR2, MAGR, MAGRP, PHI, ALPHA, BETA, COSBETA) & 297 !$OMP PRIVATE(LBRAINC, LBRA, MBRA, L, LKETINC, LKET, MKET) & 298 !$OMP PRIVATE(BLOCH, KDOTL) & 299 !$OMP PRIVATE(RCUTTB, IBRA, IKET, AMMBRA, WIGLBRAMBRA, ANGFACTOR, MP) 300 301 DO I = 1, NATS 302 303 ! Build the lists of orbitals on each atom 304 305 SELECT CASE(BASIS(ELEMPOINTER(I))) 306 307 CASE("s") 308 BASISI(1) = 0 309 BASISI(2) = -1 310 CASE("p") 311 BASISI(1) = 1 312 BASISI(2) = -1 313 CASE("d") 314 BASISI(1) = 2 315 BASISI(2) = -1 316 CASE("f") 317 BASISI(1) = 3 318 BASISI(2) = -1 319 CASE("sp") 320 BASISI(1) = 0 321 BASISI(2) = 1 322 BASISI(3) = -1 323 CASE("sd") 324 BASISI(1) = 0 325 BASISI(2) = 2 326 BASISI(3) = -1 327 CASE("sf") 328 BASISI(1) = 0 329 BASISI(2) = 3 330 BASISI(3) = -1 331 CASE("pd") 332 BASISI(1) = 1 333 BASISI(2) = 2 334 BASISI(3) = -1 335 CASE("pf") 336 BASISI(1) = 1 337 BASISI(2) = 3 338 BASISI(3) = -1 339 CASE("df") 340 BASISI(1) = 2 341 BASISI(2) = 3 342 BASISI(3) = -1 343 CASE("spd") 344 BASISI(1) = 0 345 BASISI(2) = 1 346 BASISI(3) = 2 347 BASISI(4) = -1 348 CASE("spf") 349 BASISI(1) = 0 350 BASISI(2) = 1 351 BASISI(3) = 3 352 BASISI(4) = -1 353 CASE("sdf") 354 BASISI(1) = 0 355 BASISI(2) = 2 356 BASISI(3) = 3 357 BASISI(4) = -1 358 CASE("pdf") 359 BASISI(1) = 1 360 BASISI(2) = 2 361 BASISI(3) = 3 362 BASISI(4) = -1 363 CASE("spdf") 364 BASISI(1) = 0 365 BASISI(2) = 1 366 BASISI(3) = 2 367 BASISI(4) = 3 368 BASISI(5) = -1 369 END SELECT 370 371 INDI = MATINDLIST(I) 372 373 ! open loop over neighbors J of atom I 374 DO NEWJ = 1, TOTNEBTB(I) 375 376 J = NEBTB(1, NEWJ, I) 377 378 ! IF ( J .GE. I ) THEN 379 380 PBCI = NEBTB(2, NEWJ, I) 381 PBCJ = NEBTB(3, NEWJ, I) 382 PBCK = NEBTB(4, NEWJ, I) 383 384 RIJ(1) = CR(1,J) + REAL(PBCI) * BOXDIMS(1) & 385 - CR(1,I) 386 RIJ(2) = CR(2,J) + REAL(PBCJ) * BOXDIMS(2) & 387 - CR(2,I) 388 RIJ(3) = CR(3,J) + REAL(PBCK) * BOXDIMS(3) & 389 - CR(3,I) 390 391 MAGR2 = RIJ(1)*RIJ(1) + RIJ(2)*RIJ(2) + RIJ(3)*RIJ(3) 392 393 RCUTTB = ZERO 394 395 DO K = 1, NOINT 396 397 IF ( (ATELE(I) .EQ. ELE1(K) .AND. & 398 ATELE(J) .EQ. ELE2(K)) .OR. & 399 (ATELE(J) .EQ. ELE1(K) .AND. & 400 ATELE(I) .EQ. ELE2(K) )) THEN 401 402 IF (BOND(8,K) .GT. RCUTTB ) RCUTTB = BOND(8,K) 403 404 IF (BASISTYPE .EQ. "NONORTHO") THEN 405 IF (OVERL(8,K) .GT. RCUTTB ) RCUTTB = OVERL(8,K) 406 ENDIF 407 408 ENDIF 409 410 ENDDO 411 412 IF (MAGR2 .LT. RCUTTB*RCUTTB) THEN 413 414 ! KDOTL = KPOINT(1)*RIJ(1) + KPOINT(2)*RIJ(2) + & 415 ! KPOINT(3)*RIJ(3) 416 417 ! BLOCH = EXP(CMPLX(ZERO,KDOTL)) 418 419 MAGR = SQRT(MAGR2) 420 421 SELECT CASE(BASIS(ELEMPOINTER(J))) 422 CASE("s") 423 BASISJ(1) = 0 424 BASISJ(2) = -1 425 CASE("p") 426 BASISJ(1) = 1 427 BASISJ(2) = -1 428 CASE("d") 429 BASISJ(1) = 2 430 BASISJ(2) = -1 431 CASE("f") 432 BASISJ(1) = 3 433 BASISJ(2) = -1 434 CASE("sp") 435 BASISJ(1) = 0 436 BASISJ(2) = 1 437 BASISJ(3) = -1 438 CASE("sd") 439 BASISJ(1) = 0 440 BASISJ(2) = 2 441 BASISJ(3) = -1 442 CASE("sf") 443 BASISJ(1) = 0 444 BASISJ(2) = 3 445 BASISJ(3) = -1 446 CASE("pd") 447 BASISJ(1) = 1 448 BASISJ(2) = 2 449 BASISJ(3) = -1 450 CASE("pf") 451 BASISJ(1) = 1 452 BASISJ(2) = 3 453 BASISJ(3) = -1 454 CASE("df") 455 BASISJ(1) = 2 456 BASISJ(2) = 3 457 BASISJ(3) = -1 458 CASE("spd") 459 BASISJ(1) = 0 460 BASISJ(2) = 1 461 BASISJ(3) = 2 462 BASISJ(4) = -1 463 CASE("spf") 464 BASISJ(1) = 0 465 BASISJ(2) = 1 466 BASISJ(3) = 3 467 BASISJ(4) = -1 468 CASE("sdf") 469 BASISJ(1) = 0 470 BASISJ(2) = 2 471 BASISJ(3) = 3 472 BASISJ(4) = -1 473 CASE("pdf") 474 BASISJ(1) = 1 475 BASISJ(2) = 2 476 BASISJ(3) = 3 477 BASISJ(4) = -1 478 CASE("spdf") 479 BASISJ(1) = 0 480 BASISJ(2) = 1 481 BASISJ(3) = 2 482 BASISJ(4) = 3 483 BASISJ(5) = -1 484 END SELECT 485 486 INDJ = MATINDLIST(J) 487 488 MAGRP = SQRT(RIJ(1) * RIJ(1) + RIJ(2) * RIJ(2)) 489 490 ! transform to system in which z-axis is aligned with RIJ, 491 IF (ABS(RIJ(1)) .GT. 1E-12) THEN 492 493 IF (RIJ(1) .GT. ZERO .AND. RIJ(2) .GE. ZERO) THEN 494 PHI = ZERO 495 ELSEIF (RIJ(1) .GT. ZERO .AND. RIJ(2) .LT. ZERO) THEN 496 PHI = TWO * PI 497 ELSE 498 PHI = PI 499 ENDIF 500 ALPHA = ATAN(RIJ(2) / RIJ(1)) + PHI 501 502 ELSEIF (ABS(RIJ(2)) .GT. 1E-12) THEN 503 504 IF (RIJ(2) .GT. 1E-12) THEN 505 ALPHA = PI / TWO 506 ELSE 507 ALPHA = THREE * PI / TWO 508 ENDIF 509 510 ELSE 511 ! pathological case: beta=0 and alpha undefined, but 512 ! this doesn't matter for matrix elements 513 514 ALPHA = ZERO 515 516 ENDIF 517 518 COSBETA = RIJ(3)/MAGR 519 BETA = ACOS(RIJ(3) / MAGR) 520 521 ! Build matrix elements using eqns (1)-(9) in PRB 72 165107 522 523 ! The loops over LBRA and LKET need to take into account 524 ! the orbitals assigned to each atom, e.g., sd rather than 525 ! spd... 526 527 IBRA = INDI + 1 528 529 LBRAINC = 1 530 DO WHILE (BASISI(LBRAINC) .NE. -1) 531 532 LBRA = BASISI(LBRAINC) 533 LBRAINC = LBRAINC + 1 534 535 DO MBRA = -LBRA, LBRA 536 537 ! We can calculate these two outside the 538 ! MKET loop... 539 540 AMMBRA = AM(MBRA, ALPHA) 541 WIGLBRAMBRA = WIGNERD(LBRA, ABS(MBRA), 0, COSBETA) 542 543 IKET = INDJ + 1 544 545 LKETINC = 1 546 DO WHILE (BASISJ(LKETINC) .NE. -1) 547 548 LKET = BASISJ(LKETINC) 549 LKETINC = LKETINC + 1 550 551 DO MKET = -LKET, LKET 552 553 ! This is the sigma bonds (mp = 0) 554 555 ! Hamiltonian build 556 557 ! Pre-compute the angular part so we can use it 558 ! again later if we're building the S matrix too 559 560 ANGFACTOR = TWO * AMMBRA * & 561 AM(MKET, ALPHA) * & 562 WIGLBRAMBRA * & 563 WIGNERD(LKET, ABS(MKET), 0, COSBETA) 564 565 KHTMP = CMPLX(ANGFACTOR * & 566 UNIVSCALE(I, J, LBRA, LKET, & 567 0, MAGR, "H")) 568 569 KCOUNT = 0 570 571 DO KX = 1, NKX 572 KPOINT(1) = KX0 + TWO*PI*REAL(KX-1)/(BOXDIMS(1)*REAL(NKX)) 573 DO KY = 1, NKY 574 KPOINT(2) = KY0 + TWO*PI*REAL(KY-1)/(BOXDIMS(2)*REAL(NKY)) 575 DO KZ = 1, NKZ 576 KPOINT(3) = KZ0 + TWO*PI*REAL(KZ-1)/(BOXDIMS(3)*REAL(NKZ)) 577 578 KCOUNT = KCOUNT+1 579 580 KDOTL = KPOINT(1)*RIJ(1) + KPOINT(2)*RIJ(2) + & 581 KPOINT(3)*RIJ(3) 582 583 BLOCH = EXP(CMPLX(ZERO,KDOTL)) 584 585 HK(IBRA, IKET, KCOUNT) = & 586 HK(IBRA, IKET, KCOUNT) + & 587 BLOCH*KHTMP 588 589 ENDDO 590 ENDDO 591 ENDDO 592 593 594 ! Overlap matrix build 595 596 ! IF (BASISTYPE .EQ. "NONORTHO") THEN 597 598 ! SMAT(IBRA, IKET) = SMAT(IBRA, IKET) + & 599 ! BLOCH*ANGFACTOR * & 600 ! UNIVSCALE(I, J, LBRA, LKET, 0, MAGR, "S") 601 602 ! ENDIF 603 604 ! everything else 605 606 DO MP = 1, MIN(LBRA, LKET) 607 608 ANGFACTOR = & 609 SLMMP(LBRA, MBRA, MP, ALPHA, COSBETA)* & 610 SLMMP(LKET, MKET, MP, ALPHA, COSBETA)+ & 611 TLMMP(LBRA, MBRA, MP, ALPHA, COSBETA)* & 612 TLMMP(LKET, MKET, MP, ALPHA, COSBETA) 613 614 KHTMP = CMPLX(ANGFACTOR * & 615 UNIVSCALE(I, J, LBRA, LKET, & 616 MP, MAGR, "H")) 617 618 KCOUNT = 0 619 DO KX = 1, NKX 620 KPOINT(1) = KX0 + TWO*PI*REAL(KX-1)/(BOXDIMS(1)*REAL(NKX)) 621 DO KY = 1, NKY 622 KPOINT(2) = KY0 + TWO*PI*REAL(KY-1)/(BOXDIMS(2)*REAL(NKY)) 623 DO KZ = 1, NKZ 624 KPOINT(3) = KZ0 + TWO*PI*REAL(KZ-1)/(BOXDIMS(3)*REAL(NKZ)) 625 626 KCOUNT = KCOUNT+1 627 628 KDOTL = KPOINT(1)*RIJ(1) + KPOINT(2)*RIJ(2) + & 629 KPOINT(3)*RIJ(3) 630 631 BLOCH = EXP(CMPLX(ZERO,KDOTL)) 632 633 HK(IBRA, IKET, KCOUNT) = & 634 HK(IBRA, IKET, KCOUNT) + & 635 BLOCH*KHTMP 636 637 ENDDO 638 ENDDO 639 ENDDO 640 641 642 643 ! HK(IBRA, IKET, KCOUNT) = & 644 ! HK(IBRA, IKET, KCOUNT) + & 645 ! BLOCH*CMPLX(ANGFACTOR * & 646 ! UNIVSCALE(I, J, LBRA, LKET, & 647 ! MP, MAGR, "H")) 648 649 ! IF (BASISTYPE .EQ. "NONORTHO") THEN 650 651 ! SMAT(IBRA, IKET) = SMAT(IBRA, IKET) + & 652 ! BLOCH*ANGFACTOR * & 653 ! UNIVSCALE(I, J, LBRA, LKET, MP, MAGR, "S") 654 655 ! ENDIF 656 657 ENDDO 658 659 ! HK(IKET, IBRA, KCOUNT) = CONJG(HK(IBRA, IKET, KCOUNT)) 660 ! IF (BASISTYPE .EQ. "NONORTHO") & 661 ! SMAT(IKET, IBRA) = SMAT(IBRA, IKET) 662 663 IKET = IKET + 1 664 665 ENDDO 666 667 ENDDO 668 669 IBRA = IBRA + 1 670 671 ENDDO 672 ENDDO 673 ENDIF 674 ! ENDIF 675 ENDDO 676 ! INDI = INDI + NORBI 677 678 ENDDO 679 !$OMP END PARALLEL DO 680 ! ENDDO 681 ! ENDDO 682 ! ENDDO 683 684 ! Save the diagonal elements: it will help a lot when we add in the partial charges 685 686 DO I = 1, NKTOT 687 DO J = 1,HDIM 688 HKDIAG(J,I) = HK(J,J,I) 689 ENDDO 690 ENDDO 691 692 693 ! IF (BASISTYPE .EQ. "NONORTHO") THEN 694 695 ! H0 = H 696 ! CALL GENX 697 698 ! IF (DEBUGON .EQ. 1) THEN 699 700 ! OPEN(UNIT=30, STATUS="UNKNOWN", FILE="myS.dat") 701 ! OPEN(UNIT=31, STATUS="UNKNOWN", FILE="myH0.dat") 702 703 ! PRINT*, "Caution - the Slater-Koster H and overlap matrices are being written to file" 704 705 ! DO I = 1, HDIM 706 ! WRITE(30,10) (SMAT(I,J), J = 1, HDIM) 707 ! WRITE(31,10) (H0(I,J), J = 1, HDIM) 708 ! ENDDO 709 710 ! CLOSE(30) 711 ! CLOSE(31) 712 713 !10 FORMAT(100F12.6) 714 715 ! ENDIF 716 717 ! ENDIF 718 719 RETURN 720 721END SUBROUTINE KBLDNEWH 722