1* 2* $Id$ 3* 4*********************************************************************** 5* band_minimizer * 6* * 7* This is a developing band structure parallel code for NWCHEM * 8* + tcgmsg message passing library used * 9* + my own slap-decomposed parallel 3d-FFT(real->complex) used * 10* * 11* * 12*********************************************************************** 13 14 logical function band_minimizer(rtdb,flag) 15 implicit none 16 integer rtdb 17 integer flag 18 19#include "global.fh" 20#include "bafdecls.fh" 21#include "inp.fh" 22#include "btdb.fh" 23#include "stdio.fh" 24#include "util.fh" 25#include "errquit.fh" 26 27* **** parallel variables **** 28 integer taskid,np,np_i,np_j,np_k 29 integer MASTER 30 parameter(MASTER=0) 31 32* **** timing variables **** 33 real*8 cpu1,cpu2,cpu3,cpu4 34 real*8 t1,t2,t3,t4,av 35 36* **** lattice variables **** 37 integer ngrid(3),nwave,nfft3d 38 real*8 a,b,c,alpha,beta,gamma 39 40* ***** energy variables **** 41 logical spin_orbit 42 integer ispin,ne(2) 43 real*8 E(10) 44 real*8 dipole(3) 45 real*8 stress(3,3),lstress(6) 46 47* **** gradient variables **** 48 integer fion(2) 49 50* **** error variables **** 51 logical value 52 integer ierr 53 54* **** local variables **** 55 logical newpsi,lprint,mprint,hprint 56 real*8 gx,gy,gz,cx,cy,cz 57 real*8 EV,pi,weight 58 real*8 icharge,en(2) 59 real*8 f0,f1,f2,f3,f4,f5,f6 60 integer if1,if2 61 integer i,k,ia,nion,vers,nbrillioun,nb,ms 62 integer mapping,minimizer 63 64* **** external functions **** 65 real*8 lattice_omega,lattice_unita,lattice_ecut,lattice_wcut 66 real*8 lattice_unitg,ion_amass,ion_TotalCharge 67 character*4 ion_aname 68 external lattice_omega,lattice_unita,lattice_ecut,lattice_wcut 69 external lattice_unitg,ion_amass,ion_TotalCharge 70 external ion_aname 71 72 73 real*8 control_tole,control_tolc,control_tolr,ion_rion 74 external control_tole,control_tolc,control_tolr,ion_rion 75 real*8 control_time_step,control_fake_mass 76 external control_time_step,control_fake_mass 77 logical control_read,control_move,ion_init,ion_q_FixIon 78 external control_read,control_move,ion_init,ion_q_FixIon 79 logical control_spin_orbit 80 external control_spin_orbit 81 integer control_it_in,control_it_out,control_gga,control_version 82 integer control_ngrid,pack_nwave,ion_nion,ion_natm,ion_katm 83 integer ion_nkatm 84 external control_it_in,control_it_out,control_gga,control_version 85 external control_ngrid,pack_nwave,ion_nion,ion_natm,ion_katm 86 external ion_nkatm 87 88 character*12 control_boundry 89 external control_boundry 90 91 logical brillioun_print,control_Mulliken 92 integer brillioun_nbrillioun,brillioun_nbrillq 93 integer brillioun_nbrillioun0 94 real*8 brillioun_weight_brdcst,brillioun_ks_brdcst 95 real*8 brillioun_k_brdcst,brillioun_weight 96 external brillioun_print,control_Mulliken 97 external brillioun_nbrillioun,brillioun_nbrillq 98 external brillioun_nbrillioun0 99 external brillioun_weight_brdcst,brillioun_ks_brdcst 100 external brillioun_k_brdcst,brillioun_weight 101 integer c_electron_count,linesearch_count 102 external c_electron_count,linesearch_count 103 104 real*8 nwpw_timing 105 external nwpw_timing 106 integer Cram_nwave_all_brdcst,Cram_nwave_brdcst 107 external Cram_nwave_all_brdcst,Cram_nwave_brdcst 108 109 integer ewald_ncut 110 real*8 ewald_rcut,ewald_mandelung,ewald_e 111 external ewald_ncut 112 external ewald_rcut,ewald_mandelung,ewald_e 113 logical cpsp_semicore,psi_filefind,cpsi_initialize,cpsi_finalize 114 external cpsp_semicore,psi_filefind,cpsi_initialize,cpsi_finalize 115 real*8 c_cgsd_energy,c_cgsd_noit_energy 116 external c_cgsd_energy,c_cgsd_noit_energy 117 integer cpsp_lmax,cpsp_locp,cpsp_lmmax 118 external cpsp_lmax,cpsp_locp,cpsp_lmmax 119 integer cpsp_nprj,cpsp_psp_type 120 external cpsp_nprj,cpsp_psp_type 121 real*8 cpsp_rcore,cpsp_rc,cpsp_ncore,cpsp_zv 122 external cpsp_rcore,cpsp_rc,cpsp_ncore,cpsp_zv 123 character*4 ion_atom 124 external ion_atom 125 126 character*255 cpsp_comment,comment 127 external cpsp_comment 128 129 integer cpsi_ispin,cpsi_ne,psi_get_version 130 external cpsi_ispin,cpsi_ne,psi_get_version 131 logical pspw_reformat_c_wvfnc,control_check_charge_multiplicity 132 external pspw_reformat_c_wvfnc,control_check_charge_multiplicity 133 integer control_mapping,control_minimizer,control_np_dimensions 134 external control_mapping,control_minimizer,control_np_dimensions 135 integer control_scf_algorithm 136 external control_scf_algorithm 137 real*8 control_ks_alpha,control_kerker_g0 138 external control_ks_alpha,control_kerker_g0 139 integer control_ks_maxit_orb,control_ks_maxit_orbs 140 external control_ks_maxit_orb,control_ks_maxit_orbs 141 real*8 cpsi_occupation,control_fractional_temperature 142 external cpsi_occupation,control_fractional_temperature 143 logical control_fractional,control_print 144 external control_fractional,control_print 145 integer control_fractional_smeartype 146 external control_fractional_smeartype 147 real*8 control_fractional_kT,control_fractional_alpha 148 external control_fractional_kT,control_fractional_alpha 149 150c character*255 cpsp_comment,comment 151c external cpsp_comment 152 153 integer ion_nconstraints,ion_ndof 154 external ion_nconstraints,ion_ndof 155 logical control_smooth_cutoff 156 external control_smooth_cutoff 157 real*8 control_smooth_cutoff_values 158 external control_smooth_cutoff_values 159 160 161*****************************| PROLOGUE |**************************** 162 163 value = .true. 164 pi = 4.0d0*datan(1.0d0) 165 166 call nwpw_timing_init() 167 call dcopy(10,0.0d0,0,E,1) 168 169 170* **** get parallel variables **** 171 call Parallel_Init() 172 call Parallel_np(np) 173 call Parallel_taskid(taskid) 174 175 value = control_read(5,rtdb) 176 if (.not. value) 177 > call errquit('error reading control',0, DISK_ERR) 178 179 lprint = ((taskid.eq.MASTER).and.(control_print(print_low))) 180 mprint = ((taskid.eq.MASTER).and.(control_print(print_medium))) 181 hprint = ((taskid.eq.MASTER).and.(control_print(print_high))) 182 183 if (taskid.eq.MASTER) call current_second(cpu1) 184* ***** print out header **** 185 if (mprint) then 186 write(luout,1000) 187 write(luout,1010) 188 write(luout,1020) 189 write(luout,1010) 190 write(luout,1030) 191 write(luout,1010) 192 write(luout,1035) 193 write(luout,1010) 194 write(luout,1040) 195 write(luout,1010) 196 write(luout,1041) 197 write(luout,1042) 198 write(luout,1043) 199 write(luout,1044) 200 write(luout,1010) 201 write(luout,1000) 202 call nwpw_message(1) 203 write(luout,1110) 204 call flush(luout) 205 end if 206 207 208 call Parallel3d_Init(control_np_dimensions(2), 209 > control_np_dimensions(3)) 210 call Parallel3d_np_i(np_i) 211 call Parallel3d_np_j(np_j) 212 call Parallel3d_np_k(np_k) 213 214 ngrid(1) = control_ngrid(1) 215 ngrid(2) = control_ngrid(2) 216 ngrid(3) = control_ngrid(3) 217 nwave = 0 218 minimizer = control_minimizer() 219 mapping = control_mapping() 220 ierr = 0 221 222* **** initialize C3dB data structure **** 223 call C3dB_Init(1,ngrid(1),ngrid(2),ngrid(3),mapping) 224 call C3dB_nfft3d(1,nfft3d) 225 226 call cpsi_data_init(20) 227 228* **** read ions **** 229 value = ion_init(rtdb) 230 call center_geom(cx,cy,cz) 231 call center_mass(gx,gy,gz) 232 233* **** initialize lattice data structure **** 234 call lattice_init() 235 call c_G_init() 236 237* **** nitialize brillion *** 238 call brillioun_init() 239 call Cram_Init() 240 call C3dB_pfft_init() 241 242* **** initialize D3dB data structure and mask for GGA **** 243 if ((control_gga().ge.10).and.(control_gga().le.200)) THEN 244 call D3dB_Init(1,ngrid(1),ngrid(2),ngrid(3),mapping) 245 call G_init() 246 call mask_init() 247 end if 248 249c* **** read ions **** 250c value = ion_init(rtdb) 251c call center_geom(cx,cy,cz) 252c call center_mass(gx,gy,gz) 253 254* **** allocate psp data structure and read in psedupotentials into it **** 255 call cpsp_init() 256 call cpsp_readall() 257 if (cpsp_semicore(0)) call c_semicore_check() 258 259 260* **** initialize ke,and coulomb data structures **** 261 call cstrfac_init() 262 call cke_init() 263 call c_coulomb_init() 264 call ewald_init() 265 266* **** set up phase factors at the current geometry **** 267 call cphafac() 268 call cphafac_k() 269 call ewald_phafac() 270 271* **** read in wavefunctions and initialize psi **** 272 if (.not.control_check_charge_multiplicity()) then 273 call cpsi_new() 274 newpsi = .true. 275 else 276 newpsi = .false. 277* **** convert from pspw format to band format **** 278 vers = psi_get_version() 279 if ((vers.eq.3).or.(vers.eq.4)) then 280 newpsi = .true. 281 value = btdb_parallel(.false.) 282 if (taskid.eq.MASTER) then 283 value= pspw_reformat_c_wvfnc(1) 284 end if 285 value = btdb_parallel(.true.) 286 end if 287 end if 288 289 call psi_get_ne(ispin,ne) 290 if (ispin.eq.3) then 291 spin_orbit = .true. 292 ispin=2 293 else 294 spin_orbit = .false. 295 end if 296 nbrillioun = brillioun_nbrillioun() 297 call Pneb_init(ispin,ne,nbrillioun,spin_orbit) 298 value = cpsi_initialize(.true.) 299 newpsi = newpsi.or.value 300 301* **** electron and geodesic data structures **** 302 call c_electron_init() 303 call c_geodesic_init() 304 305 306* **** initialize QM/MM **** 307 call band_init_APC(rtdb) 308 309 310* **** initialize HFX **** 311 call band_init_HFX(rtdb,nbrillioun,ispin,ne) 312 313 314* **** initialize FixIon constraint **** 315 call ion_init_FixIon(rtdb) 316 317* **** initialize linesearching **** 318 call linesearch_init() 319 320 321* |**************************| 322****************** summary of input data ********************** 323* |**************************| 324 325 326* **** determine en **** 327 if (.not.control_spin_orbit()) then 328 icharge = 0.0d0 329 en(1) = 0.0d0 330 en(2) = 0.0d0 331 b = dble(3-cpsi_ispin()) 332 do nb=1,brillioun_nbrillq() 333 weight = brillioun_weight(nb) 334 do ms=1,cpsi_ispin() 335 do i=1,ne(ms) 336 a = cpsi_occupation(nb,ms,i) 337 icharge = icharge - b*a*weight 338 en(ms) = en(ms) + a*weight 339 end do 340 end do 341 end do 342 call K1dB_Vector_SumAll(2,en) 343 call K1dB_SumAll(icharge) 344 else 345 icharge = -cpsi_ne(1) 346 en(1) = cpsi_ne(1) 347 en(cpsi_ispin()) = cpsi_ne(cpsi_ispin()) 348 end if 349 350 351 if (mprint) then 352 write(luout,1111) np 353 write(luout,1117) np_i,np_j,np_k 354 if (mapping.eq.1) write(luout,1112) 355 if (mapping.eq.2) write(luout,1113) 356 if (mapping.eq.3) write(luout,1118) 357 358 write(luout,1115) 359 write(luout,1121) control_boundry(),control_version() 360 if (.not.spin_orbit) then 361 if (cpsi_ispin().eq.1) write(luout,1130) "restricted" 362 if (cpsi_ispin().eq.2) write(luout,1130) "unrestricted" 363 else 364 write(luout,1130) "spin orbit" 365 end if 366 367 call v_bwexc_print(luout,control_gga()) 368 369 call band_print_HFX(luout) 370 write(luout,1140) 371 do ia = 1,ion_nkatm() 372 write(luout,1150) ia,ion_atom(ia), 373 > cpsp_zv(ia),cpsp_lmax(ia) 374 comment = cpsp_comment(ia) 375 i = inp_strlen(comment) 376 write(6,1157) comment(1:i) 377 write(6,1158) cpsp_psp_type(ia) 378 write(luout,1152) cpsp_lmax(ia) 379 write(luout,1153) cpsp_locp(ia) 380c write(luout,1154) cpsp_lmmax(ia) 381 write(luout,1154) cpsp_nprj(ia) 382 if (cpsp_semicore(ia)) 383 > write(luout,1155) cpsp_rcore(ia),cpsp_ncore(ia) 384 write(luout,1151) (cpsp_rc(i,ia),i=0,cpsp_lmax(ia)) 385 end do 386 387 icharge = icharge + ion_TotalCharge() 388 write(luout,1159) icharge 389 390 391 write(luout,1160) 392 write(luout,1170) (ion_atom(K),ion_natm(K),K=1,ion_nkatm()) 393 if (hprint) then 394 write(luout,1180) 395 do I=1,ion_nion() 396 if (ion_q_FixIon(I)) then 397 write(luout,1191) I,ion_aname(I),(ion_rion(K,I),K=1,3), 398 > ion_amass(I)/1822.89d0 399 else 400 write(luout,1190) I,ion_aname(I),(ion_rion(K,I),K=1,3), 401 > ion_amass(I)/1822.89d0 402 end if 403 end do 404 write(luout,1200) cx,cy,cz 405 write(luout,1210) gx,gy,gz 406 write(luout,1211) ion_nconstraints(),ion_ndof() 407 endif 408 write(luout,1220) en(1),en(cpsi_ispin()), 409 > ' (Fourier space)' 410 write(luout,1221) cpsi_ne(1),cpsi_ne(cpsi_ispin()), 411 > ' (Fourier space)' 412 413 write(luout,1230) 414 write(luout,1241) lattice_unita(1,1), 415 > lattice_unita(2,1), 416 > lattice_unita(3,1) 417 write(luout,1242) lattice_unita(1,2), 418 > lattice_unita(2,2), 419 > lattice_unita(3,2) 420 write(luout,1243) lattice_unita(1,3), 421 > lattice_unita(2,3), 422 > lattice_unita(3,3) 423 write(luout,1244) lattice_unitg(1,1), 424 > lattice_unitg(2,1), 425 > lattice_unitg(3,1) 426 write(luout,1245) lattice_unitg(1,2), 427 > lattice_unitg(2,2), 428 > lattice_unitg(3,2) 429 write(luout,1246) lattice_unitg(1,3), 430 > lattice_unitg(2,3), 431 > lattice_unitg(3,3) 432 call lattice_abc_abg(a,b,c,alpha,beta,gamma) 433 write(luout,1232) a,b,c,alpha,beta,gamma 434 write(luout,1231) lattice_omega() 435 write(luout,1260) ewald_rcut(),ewald_ncut() 436 write(luout,1261) ewald_mandelung() 437 438 write(luout,1255) 439 write(luout,1256) brillioun_nbrillioun(), 440 > brillioun_nbrillioun0() 441 end if 442 443c **** print brillioun zone - extra logic for distributed kpoints **** 444 if (brillioun_print()) then 445 do i=1,brillioun_nbrillioun() 446 f0 = brillioun_weight_brdcst(i) 447 f1 = brillioun_ks_brdcst(1,i) 448 f2 = brillioun_ks_brdcst(2,i) 449 f3 = brillioun_ks_brdcst(3,i) 450 f4 = brillioun_k_brdcst(1,i) 451 f5 = brillioun_k_brdcst(2,i) 452 f6 = brillioun_k_brdcst(3,i) 453 if (mprint) write(luout,1257) f0,f1,f2,f3,f4,f5,f6 454 end do 455 else 456 if (mprint) write(luout,1258) 457 end if 458 459 if1 = Cram_nwave_all_brdcst(0) 460 if2 = Cram_nwave_brdcst(0) 461 if (mprint) then 462 write(luout,1249) 463 write(luout,1250) lattice_ecut(),ngrid(1),ngrid(2),ngrid(3), 464 > if1,if2 465 end if 466 467 if (brillioun_print()) then 468 do i=1,brillioun_nbrillioun() 469 if1 = Cram_nwave_all_brdcst(i) 470 if2 = Cram_nwave_brdcst(i) 471 if (mprint) then 472 write(luout,1251) i,lattice_wcut(),ngrid(1),ngrid(2),ngrid(3), 473 > if1,if2 474 end if 475 end do 476 else 477 if (mprint) write(luout,1252) lattice_wcut() 478 end if 479 480 if (mprint) then 481 if (control_smooth_cutoff()) 482 > write(luout,1253) 483 > control_smooth_cutoff_values(1)*lattice_wcut(), 484 > control_smooth_cutoff_values(2) 485 486 write(luout,1270) 487 write(luout,1280) control_time_step(),control_fake_mass() 488 write(luout,1290) control_tole(),control_tolc() 489 write(luout,1281) control_it_in()*control_it_out(), 490 > control_it_in(),control_it_out() 491 492 if ((minimizer.eq.5).or.(minimizer.eq.8)) then 493 write(6,1291) 494 write(6,1292) 495 write(luout,1295) control_ks_maxit_orb(), 496 > control_ks_maxit_orbs() 497 if (control_scf_algorithm().eq.0) 498 > write(6,1293) "simple mixing" 499 if (control_scf_algorithm().eq.1) 500 > write(6,1293) "Anderson potential mixing" 501 if (control_scf_algorithm().eq.2) 502 > write(6,1293) "Johnson-Pulay mixing" 503 if (control_scf_algorithm().eq.3) 504 > write(6,1293) "Anderson density mixing" 505 if (minimizer.eq.5) write(luout,1296) "potential" 506 if (minimizer.eq.8) write(luout,1296) "density" 507 write(6,1294) control_ks_alpha() 508 write(6,1301) control_kerker_g0() 509 end if 510 if (control_fractional()) then 511 write(6,1297) 512 if (control_fractional_smeartype().eq.0) 513 > write(6,1298) "step function" 514 if (control_fractional_smeartype().eq.1) 515 > write(6,1298) "Fermi-Dirac" 516 if (control_fractional_smeartype().eq.2) 517 > write(6,1298) "Gaussian" 518 write(6,1299) control_fractional_kT(), 519 > control_fractional_temperature(), 520 > control_fractional_alpha() 521 end if 522 write(luout,1300) 523 call util_flush(luout) 524 end if 525 526 527 528* |***************************| 529****************** call CG minimizer ********************** 530* |***************************| 531 if (taskid.eq.MASTER) call current_second(cpu2) 532 533 534* **** calculate energy **** 535 if (flag.eq.-1) then 536 537 EV= c_cgsd_noit_energy() 538 else 539 540 EV= c_cgsd_energy(newpsi) 541 end if 542 543* **** calculate excited state orbitals **** 544 call ga_sync() 545 call c_cgsd_excited() 546 547* **** extra energy output for QA test **** 548 if (lprint) write(luout,1600) EV 549 550* **** calculate the spin contamination **** 551 if (flag.gt.-1) call cpsi_spin2(dipole(1)) 552 553* **** calculate the dipole *** 554 dipole(1) = 0.0d0 555 dipole(2) = 0.0d0 556 dipole(3) = 0.0d0 557 558* **** calculate gradient *** 559 if (flag.gt.0) then 560 nion = ion_nion() 561 value = BA_push_get(mt_dbl,(3*nion), 562 > 'fion',fion(2),fion(1)) 563 if (.not. value) call errquit('out of stack memory',0, MA_ERR) 564 565 566 call c_cgsd_energy_gradient(dbl_mb(fion(1))) 567 call dscal(3*nion,(-1.0d0),dbl_mb(fion(1)),1) 568 end if 569 570* **** calculate the stress tensor **** 571 call dcopy(9,0.0d0,0,stress,1) 572 call dcopy(6,0.0d0,0,lstress,1) 573 if (flag.eq.3) then 574 575 call cpsp_stress_init() 576 call cpsp_stress_readall() 577 call c_cgsd_energy_stress(stress,lstress) 578 call cpsp_stress_end() 579 end if 580 581 582* ************************************************************* 583* **** output energy, dipole, and gradient to rtdb for use **** 584* **** by task_energy and task_gradient **** 585* ************************************************************* 586 value = .true. 587 if (flag.gt.-1) then 588 value = btdb_put(rtdb,'band:energy',mt_dbl,1,EV) 589 value = value.and. 590 > btdb_put(rtdb,'band:dipole',mt_dbl, 591 > 3,dipole) 592 end if 593 if (flag.gt.0) then 594 value = value.and.btdb_put(rtdb,'band:gradient',mt_dbl, 595 > 3*nion,dbl_mb(fion(1))) 596 value = value.and.BA_pop_stack(fion(2)) 597 end if 598 if (flag.eq.3) then 599 value = value.and. 600 > btdb_put(rtdb,'band:stress',mt_dbl, 601 > 9,stress) 602 value = value.and. 603 > btdb_put(rtdb,'band:lstress',mt_dbl, 604 > 6,lstress) 605 end if 606 if (.not. value) 607 > call errquit('band_minimizer: error writing rtdb',0, RTDB_ERR) 608 609 610 611 if (taskid.eq.MASTER) call current_second(cpu3) 612 613* |***************************| 614****************** Epilogue ********************** 615* |***************************| 616 617* **** calculate interpolated band structure plot *** 618 call band_interpolate_structure(rtdb) 619cc call band_kp_structure(rtdb) 620 621* **** calculate Mulliken Populations *** 622 if (control_Mulliken()) then 623 call cpsi_Mulliken(rtdb) 624 end if 625 626 627* **** write geometry to rtdb **** 628 call ion_write(rtdb) 629 630* **** write wavefunctions to file and finalize psi **** 631 if (flag.eq.-1) then 632 value = cpsi_finalize(.false.) 633 else 634 value = cpsi_finalize(.true.) 635 end if 636 637* **** deallocate heap memory **** 638 call c_electron_finalize() 639 call c_geodesic_finalize() 640 call ewald_end() 641 call cstrfac_end() 642 call c_coulomb_end() 643 call cke_end() 644 call cpsp_end() 645 call Cram_end() 646 call c_G_end() 647 call band_end_HFX() 648 call brillioun_end() 649 call band_end_APC() 650 call ion_end() 651 call ion_end_FixIon() 652 call cpsi_data_end() 653 call C3dB_pfft_end() 654 call C3dB_end(1) 655 IF ((control_gga().ge.10).and.(control_gga().le.200)) THEN 656 call mask_end() 657 call G_end() 658 call D3dB_end(1) 659 end if 660 661* |***************************| 662****************** report consumed cputime ********************** 663* |***************************| 664 if (lprint) then 665 CALL current_second(cpu4) 666 667 T1=CPU2-CPU1 668 T2=CPU3-CPU2 669 T3=CPU4-CPU3 670 T4=CPU4-CPU1 671 AV=T2/dble(c_electron_count()) 672 write(luout,1801) 673 write(luout,1802) 674 write(luout,1803) T1 675 write(luout,1804) T2 676 write(luout,1805) T3 677 write(luout,1806) T4 678 write(luout,1807) AV,c_electron_count(),linesearch_count() 679 680 call nwpw_timing_print_final(mprint,c_electron_count()) 681 write(luout,*) 682 CALL nwpw_MESSAGE(4) 683 end if 684 685 686 call Parallel3d_Finalize() 687 call Parallel_Finalize() 688 band_minimizer = value 689 return 690 691 692*::::::::::::::::::::::::::: format ::::::::::::::::::::::::::::::::: 693 1000 FORMAT(10X, 694 > '**********************************************************') 695 1010 FORMAT(10X, 696 > '* *') 697 1020 FORMAT(10X, 698 > '* NWPW BAND Calculation *') 699 1030 FORMAT(10X, 700 > '* [(bundled Grassmann/Stiefel manifold implementation)] *') 701 1035 FORMAT(10x, 702 > '* [ NorthWest Chemistry implementation ] *') 703 1040 FORMAT(10X, 704 > '* version #1.10 01/31/03 *') 705 1041 FORMAT(10X, 706 > '* A pseudopotential plane-wave band structure program *') 707 1042 FORMAT(10X, 708 > '* with Brillouin zone sampling for optimizing crystals, *') 709 1043 FORMAT(10X, 710 > '* slabs, and polymers. Developed by Eric J. Bylaska, *') 711 1044 FORMAT(10X, 712 > '* Edoardo Apra, and Patrick Nichols. *') 713 1100 FORMAT(//) 714 1110 FORMAT(10X,'================ input data ========================') 715 1111 FORMAT(/' number of processors used:',I16) 716 1112 FORMAT( ' parallel mapping : 1d slab') 717 1113 FORMAT( ' parallel mapping : 2d hilbert') 718 1115 FORMAT(/' options:') 719 1117 FORMAT( ' processor grid :',I4,' x',I4,' x',I4) 720 1118 FORMAT( ' parallel mapping : 2d hcurve') 721 1120 FORMAT(5X,' ionic motion = ',A) 722 1121 FORMAT(5X,' boundary conditions = ',A,'(version', I1,')') 723 1130 FORMAT(5X,' electron spin = ',A) 724 1131 FORMAT(5X,' exchange-correlation = ',A) 725 1140 FORMAT(/' elements involved in the cluster:') 726 1150 FORMAT(5X,I2,': ',A4,' core charge:',F4.1,' lmax=',I1) 727 1151 FORMAT(5X,' cutoff =',4F8.3) 728 1152 FORMAT(12X,' highest angular component : ',i3) 729 1153 FORMAT(12X,' local potential used : ',i3) 730 1154 FORMAT(12X,' number of non-local projections: ',i3) 731 1155 FORMAT(12X,' semicore corrections included : ', 732 > F6.3,' (radius) ',F6.3,' (charge)') 733 1156 FORMAT(12X,' aperiodic cutoff radius : ',F6.3) 734 1157 FORMAT(12X,' comment : ',A) 735 1158 FORMAT(12X,' pseudpotential type : ',i3) 736 737 1159 FORMAT(/' total charge:',F8.3) 738 1160 FORMAT(/' atomic composition:') 739 1170 FORMAT(7(5X,A4,':',I3)) 740 1180 FORMAT(/' initial position of ions:') 741 1190 FORMAT(5X, I4, A5 ,' (',3F11.5,' ) - atomic mass= ',F7.3,' ') 742 1191 FORMAT(5X, I4, A5, ' (',3F11.5, 743 > ' ) - atomic mass= ',F7.3,' - fixed') 744 1200 FORMAT(5X,' G.C. ',' (',3F11.5,' )') 745 1210 FORMAT(5X,' C.O.M.',' (',3F11.5,' )') 746 1211 FORMAT(5X,' number of constraints = ', I6,' ( DOF = ',I6,' )' ) 747 1220 FORMAT(/' number of electrons: spin up=',F8.2, 748 > ' spin down=',F8.2,A) 749 1221 FORMAT( ' number of orbitals: spin up=',I8, 750 > ' spin down=',I8,A) 751 1230 FORMAT(/' supercell:') 752 1231 FORMAT(5x,' volume : ',F10.1) 753 1232 FORMAT(/5x,' lattice: a=',f8.3,' b=',f8.3,' c=',f8.3, 754 > /5x,' alpha=',f8.3,' beta=',f8.3,' gamma=',f8.3) 755 1241 FORMAT(5x,' lattice: a1=<',3f8.3,' >') 756 1242 FORMAT(5x,' a2=<',3f8.3,' >') 757 1243 FORMAT(5x,' a3=<',3f8.3,' >') 758 1244 FORMAT(5x,' reciprocal: b1=<',3f8.3,' >') 759 1245 FORMAT(5x,' b2=<',3f8.3,' >') 760 1246 FORMAT(5x,' b3=<',3f8.3,' >') 761 762 1249 FORMAT(/' computational grids:') 763 1250 FORMAT(5X,' density cutoff=',F7.3,' fft=',I4,'x',I4,'x',I4, 764 & '( ',I8,' waves ',I8,' per task)') 765 1251 FORMAT(5X,' wavefnc ',I3,' cutoff=',F7.3, 766 & ' fft=',I4,'x',I4,'x',I4, 767 & '( ',I8,' waves ',I8,' per task)') 768 1252 FORMAT(5x,' wavefnc cutoff=',F7.3, 769 > ' wavefunction grids not printed - ', 770 > 'number of k-points is very large') 771 1253 FORMAT(5X,' smooth wavefnc cutoff: A=',F7.3,' sigma = ',F7.3) 772 773 774 1255 FORMAT(/' brillouin zone:') 775 1256 FORMAT(5x,' number of zone points:',I6, 776 > ' (reduced from ',I0,' zone points without symmetry)') 777 1257 FORMAT(5x,' weight=',f8.3,' ks=<',3f8.3,' >, k=<',3f8.3,'>') 778 1258 FORMAT(5x,' number of k-points is very large') 779 780 1260 FORMAT(5X,' Ewald summation: cut radius=',F8.2,' and',I3) 781 1261 FORMAT(5X,' madelung=',f14.8) 782 783 1270 FORMAT(/' technical parameters:') 784 1280 FORMAT(5X, ' time step=',F10.2,5X,'fictitious mass=',F10.1) 785 1281 FORMAT(5X, ' maximum iterations =',I10, 786 > ' ( ',I4,' inner ',I6,' outer )') 787 1290 FORMAT(5X, ' tolerance=',E8.3,' (energy)',E12.3, 788 & ' (density)') 789 1291 FORMAT(/' Kohn-Sham scf parameters:') 790 1292 FORMAT(5X, ' Kohn-Sham algorithm = conjugate gradient') 791 1293 FORMAT(5X, ' scf algorithm = ',A) 792 1294 FORMAT(5X, ' scf mixing parameter =',F7.4) 793 1295 FORMAT(5X, ' Kohn-Sham iterations = ',I3, 794 > ' (',I3,' outer)') 795 1296 FORMAT(5X, ' scf mixing type = ',A) 796 1297 FORMAT(/' fractional smearing parameters:') 797 1298 FORMAT(5X, ' smearing algorithm = ',A) 798 1299 FORMAT(5X, ' smearing parameter = ',E9.3,' (',F7.1,' K)'/, 799 > 5X, ' mixing parmameter = ',F7.4) 800 1300 FORMAT(//) 801 1301 FORMAT(5X, ' Kerker damping =',F7.4) 802 1305 FORMAT(10X,'================ iteration =========================') 803 1310 FORMAT(I8,E20.10,3E15.5) 804 1320 FORMAT(' number of electrons: spin up=',F11.5,' down=',F11.5,A) 805 1330 FORMAT(/' comparison between hamiltonian and lambda matrix') 806 1340 FORMAT(I3,2I3,' H=',E16.7,', L=',E16.7,', H-L=',E16.7) 807 1350 FORMAT(/' orthonormality') 808 1360 FORMAT(I3,2I3,E18.7) 809 1370 FORMAT(I3) 810 1380 FORMAT(' ''',a,'''',I4) 811 1390 FORMAT(I3) 812 1400 FORMAT(I3,3E18.8/3X,3E18.8) 813 1410 FORMAT(10X,'============= summary of results =================') 814 1420 FORMAT( ' final position of ions:') 815 1430 FORMAT(/' total energy :',E19.10,' (',E15.5,'/ion)') 816 1440 FORMAT( ' total orbital energy:',E19.10,' (',E15.5,'/electron)') 817 1450 FORMAT( ' hartree energy :',E19.10,' (',E15.5,'/electron)') 818 1460 FORMAT( ' exc-corr energy :',E19.10,' (',E15.5,'/electron)') 819 1470 FORMAT( ' ion-ion energy :',E19.10,' (',E15.5,'/ion)') 820 1480 FORMAT(/' K.S. kinetic energy :',E19.10,' (',E15.5,'/electron)') 821 1490 FORMAT( ' K.S. V_l energy :',E19.10,' (',E15.5,'/electron)') 822 1495 FORMAT( ' K.S. V_nl energy :',E19.10,' (',E15.5,'/electron)') 823 1496 FORMAT( ' K.S. V_Hart energy :',E19.10,' (',E15.5,'/electron)') 824 1497 FORMAT( ' K.S. V_xc energy :',E19.10,' (',E15.5,'/electron)') 825 1498 FORMAT( ' Virial Coefficient :',E19.10) 826 1500 FORMAT(/' orbital energies:') 827 1510 FORMAT(2(E18.7,' (',F8.3,'eV)')) 828 1600 FORMAT(/' Total BAND energy :',E19.10) 829 830 1801 FORMAT(//'== Timing ==') 831 1802 FORMAT(/'cputime in seconds') 832 1803 FORMAT( ' prologue : ',E14.6) 833 1804 FORMAT( ' main loop : ',E14.6) 834 1805 FORMAT( ' epilogue : ',E14.6) 835 1806 FORMAT( ' total : ',E14.6) 836 1807 FORMAT( ' cputime/step: ',E14.6, 837 > ' (',I8,' evalulations,', I8,' linesearches)') 838 1808 FORMAT(A,E14.6,E14.6) 839 1809 FORMAT(//A,2A14) 840 841 9010 FORMAT(//' >> job terminated due to code =',I3,' <<') 842 843 9000 if (taskid.eq.MASTER) write(6,9010) ierr 844 call Parallel_Finalize() 845 846 band_minimizer = value 847 return 848 END 849 850