1* 2* $Id$ 3* 4 5* *********************************** 6* * * 7* * control_read * 8* * * 9* *********************************** 10 logical function control_read(code_in,rtdb) 11 implicit none 12 integer code_in 13 integer rtdb 14 15#include "inp.fh" 16#include "bafdecls.fh" 17#include "btdb.fh" 18#include "control.fh" 19#include "nwpwxc.fh" 20#include "util.fh" 21 22 integer MASTER,taskid 23 parameter (MASTER=0) 24 25 logical value,np_default,value2,value3 26 integer ispin0,ne(2),nbrill,np1,np2 27 integer nmult(5) 28 real*8 error,thresh 29 30* **** control_nose common block **** 31 logical nose 32 real*8 Pe,Te,Pr,Tr 33 common / control_nblock / Pe,Te,Pr,Tr,nose 34 35* **** control_rtdb common block **** 36 integer trtdb 37 common / control_rtdb1 / trtdb 38 39* **** control_print common block **** 40 integer print_level 41 common / control_print1 / print_level 42 43 !character*30 cell_name 44 character*50 rtdb_unita,rtdb_unitaf,rtdb_ngrid,rtdb_boundry 45 character*50 exchange_correlation,rtdb_ngrid_small 46 integer i,l,np 47! integer ind ! unused 48 49 integer control_num_kvectors 50 external control_num_kvectors 51 logical control_allow_translation 52 external control_allow_translation 53 integer Parallel_threadid 54 external Parallel_threadid 55 56 57 value = .true. 58 59 call Parallel_taskid(taskid) 60 call nwpw_timing_start(50) 61 value2 = btdb_parallel(.true.) 62 code = code_in 63 trtdb = rtdb 64 65 call util_print_get_level(print_level) 66 !write(*,*) "Print level is ",print_level 67 68 69 DO 10, i = 1,100 70* **** get parallel mappings **** 71 if (.not.btdb_get(rtdb,'nwpw:mapping',mt_int,1,mapping)) then 72 mapping = 2 73 end if 7410 continue 75 76* **** set mapping1d **** 77 if (.not.btdb_get(rtdb,'nwpw:mapping1d',mt_int,1,mapping1d)) then 78 mapping1d = 1 79 end if 80 81* **** set np_dimensions **** 82 call Parallel_np(np) 83 84**** if NOMPI then no processor groups or communicators **** 85#ifdef NOMPI 86 np_dimensions(1) = np 87 np_dimensions(2) = 1 88 np_dimensions(3) = 1 89 np_default = .true. 90#else 91 if(.not.btdb_get(rtdb,'nwpw:np_dimensions', 92 > mt_int,3,np_dimensions)) then 93 94 np_dimensions(1) = np 95 np_dimensions(2) = 1 96 np_dimensions(3) = 1 97 np_default = .true. 98 else 99 np_default = .false. 100 if (.not.((code.eq.5) .or. 101 > (code.eq.10).or. 102 > (code.eq.13).or. 103 > (code.eq.14))) np_dimensions(3) = 1 104 if (np_dimensions(3).lt.1) np_dimensions(3) = 1 105 if (np_dimensions(2).lt.1) np_dimensions(2) = 1 106 if (np_dimensions(1).lt.1) np_dimensions(1) = 1 107 108* **** reset np_dimensions(3) if larger than nbrill **** 109 nbrill = control_num_kvectors() 110 if (np_dimensions(3).gt.nbrill) np_dimensions(3) = nbrill 111 112* **** reset np_dimensions(3) if it is not a multiple of np **** 113 do while ((mod(np,np_dimensions(3)).ne.0).and. 114 > (np_dimensions(3).gt.1)) 115 np_dimensions(3) = np_dimensions(3) - 1 116 end do 117* **** reset np_dimensions(2) if it is not a multiple of np2 **** 118 np = np/np_dimensions(3) 119 do while ((mod(np,np_dimensions(2)).ne.0).and. 120 > (np_dimensions(2).gt.1)) 121 np_dimensions(2) = np_dimensions(2) - 1 122 end do 123 124 !*** temporary restriction until ne parallelized in band *** 125 if ((code.eq.5) .or. 126 > (code.eq.10).or. 127 > (code.eq.13).or. 128 > (code.eq.14)) np_dimensions(2) = 1 129 130 np_dimensions(1) = np/(np_dimensions(2)*np_dimensions(3)) 131 132 endif 133#endif 134 135* **** get balance mapping **** 136 if (.not.btdb_get(rtdb,'nwpw:balance',mt_log,1,balance)) then 137 balance = .true. 138 end if 139 140* **** get parallel io **** 141 if (.not.btdb_get(rtdb,'nwpw:parallel_io',mt_log,1,pio)) then 142 pio = .false. 143 end if 144 !pio = pio.and.(np_dimensions(2).gt.1) 145 146* **** get fast_erf **** 147 if (.not.btdb_get(rtdb,'nwpw:fast_erf',mt_log,1,fast_erf)) then 148 fast_erf = .false. 149 end if 150 151* **** get fmm **** 152 if (.not.btdb_get(rtdb,'nwpw:fmm',mt_log,1,fmm)) then 153 fmm = .false. 154 end if 155 if (.not.btdb_get(rtdb,'nwpw:fmm_lmax',mt_int,1,fmm_lmax)) then 156 fmm_lmax = 10 157 end if 158 if (.not.btdb_get(rtdb,'nwpw:fmm_lr',mt_int,1,fmm_lr)) then 159 fmm_lr = 1 160 end if 161 162* **** get periodic_dipole **** 163 if (.not.btdb_get(rtdb,'nwpw:periodic_dipole', 164 > mt_log,1,periodic_dipole)) then 165 periodic_dipole = .false. 166 end if 167 168* **** get smooth_cuoff **** 169 smooth_cutoff = .true. 170 if (.not.btdb_get(rtdb,'nwpw:smooth_cutoff', 171 > mt_dbl,2,smooth_cutoff_values)) then 172 smooth_cutoff = .false. 173 smooth_cutoff_values(1) = 2.0d0 174 smooth_cutoff_values(2) = 4.0d0 175 end if 176 177* **** get hess_model mapping **** 178 if (.not.btdb_get(rtdb,'nwpw:hess_model',mt_log, 179 > 1,hess_model)) then 180 hess_model = .false. 181 end if 182 183 184 185* ********************************* 186* **** cpsd and band_sd: stuff **** 187* ********************************* 188 if ((code.eq.1).or.(code.eq.13)) then 189 if (.not.btdb_cget(rtdb,'cpsd:cell_name',1,cell_name)) then 190 cell_name = 'cell_default' 191 end if 192 193 if (.not.btdb_cget(rtdb,'cpsd:input_wavefunction_filename', 194 > 1,input_wavefunction_filename)) then 195 if (.not.btdb_cget(rtdb,'pspw:input vectors', 196 > 1,input_wavefunction_filename)) then 197 call util_file_prefix('movecs', 198 > input_wavefunction_filename) 199 end if 200 end if 201 202 if (.not.btdb_cget(rtdb,'cpsd:output_wavefunction_filename', 203 > 1,output_wavefunction_filename)) then 204 if (.not.btdb_cget(rtdb,'pspw:output vectors', 205 > 1,output_wavefunction_filename)) then 206 call util_file_prefix('movecs', 207 > output_wavefunction_filename) 208 end if 209 end if 210 211 212 if (.not.btdb_cget(rtdb,'cpsd:exchange_correlation', 213 > 1,exchange_correlation)) 214 > exchange_correlation = 'vosko' 215 216!$OMP single 217 if (nwpwxc_rtdb_load(rtdb,"dft")) then 218c call nwpwxc_print() 219 has_disp = nwpwxc_has_disp() 220 endif 221!$OMP end single 222 223#include "control_gga.fh" 224 225 if(.not.btdb_get(rtdb,'cpsd:geometry_optimize',mt_log,1,move)) 226 > move = .false. 227 if(.not.btdb_get(rtdb,'cpsd:fractional_coordinates', 228 > mt_log,1,frac_coord)) 229 > frac_coord = .false. 230 if(.not.btdb_get(rtdb,'cpsd:npsp',mt_int,1,npsp)) 231 > npsp=0 232 if(.not.btdb_get(rtdb,'cpsd:fake_mass',mt_dbl,1,fake_mass)) 233 > fake_mass = 400000.0d0 234 if(.not.btdb_get(rtdb,'cpsd:time_step',mt_dbl,1,time_step)) 235 > time_step = 5.8d0 236 if(.not.btdb_get(rtdb,'cpsd:loop',mt_int,2,loop)) then 237 loop(1)=10 238 loop(2)=100 239 end if 240 if(.not.btdb_get(rtdb,'cpsd:tolerances',mt_dbl,3,tolerances)) 241 >then 242 tolerances(1) = 1.0d-9 243 tolerances(2) = 1.0d-9 244 tolerances(3) = 1.0d-4 245 end if 246 scaling(1) = 0.0d0 247 scaling(2) = 0.0d0 248 249 if(.not.btdb_get(rtdb,'cpsd:ecut',mt_dbl,1,ecut)) 250 > ecut=9000.0d0 251 if(.not.btdb_get(rtdb,'cpsd:wcut',mt_dbl,1,wcut)) 252 > wcut=ecut 253 if(.not.btdb_get(rtdb,'cpsd:rcut',mt_dbl,1,rcut)) 254 > rcut = 0.0d0 255 if(.not.btdb_get(rtdb,'cpsd:ncut',mt_int,1,ncut)) 256 > ncut = 1 257 if(.not.btdb_get(rtdb,'cpsd:mult',mt_int,1,multiplicity)) 258 > multiplicity = 1 259 if(.not.btdb_get(rtdb,'cpsd:ispin',mt_int,1,ispin)) 260 > ispin=1 261 262 value3 = .false. 263 if (.not.btdb_get(rtdb,'nwpw:dof_rotation',mt_log,1,dof_rotation)) 264 > then 265 dof_rotation = .false. 266 value3= .true. 267 end if 268 269 if (.not.btdb_get(rtdb,'nwpw:rotation',mt_log,1,rotation)) then 270 rotation = .true. 271 else 272 if (value3) dof_rotation = rotation 273 end if 274 275 276 if (.not.btdb_get(rtdb,'nwpw:spin_orbit',mt_log,1,spin_orbit)) 277 > spin_orbit=.false. 278 if (spin_orbit) ispin=2 279 280 281* **************************************** 282* **** cpmd and band_cpmd code: stuff **** 283* **************************************** 284 else if ((code.eq.2) .or. (code.eq.14)) then 285 if (.not.btdb_cget(rtdb,'cpmd:cell_name',1,cell_name)) then 286 cell_name = 'cell_default' 287 end if 288 289 if(.not.btdb_cget(rtdb,'cpmd:input_wavefunction_filename', 290 > 1,input_wavefunction_filename)) then 291 if(.not.btdb_cget(rtdb,'pspw:input vectors', 292 > 1,input_wavefunction_filename)) then 293 call util_file_prefix('movecs', 294 > input_wavefunction_filename) 295 end if 296 end if 297 298 if(.not.btdb_cget(rtdb,'cpmd:output_wavefunction_filename', 299 > 1,output_wavefunction_filename)) then 300 if(.not.btdb_cget(rtdb,'pspw:output vectors', 301 > 1,output_wavefunction_filename)) then 302 call util_file_prefix('movecs', 303 > output_wavefunction_filename) 304 end if 305 end if 306 307 if(.not.btdb_cget(rtdb,'cpmd:input_v_wavefunction_filename', 308 > 1,input_v_wavefunction_filename)) then 309 if(.not.btdb_cget(rtdb,'pspw:input vvectors', 310 > 1,input_v_wavefunction_filename)) then 311 call util_file_prefix('vmovecs', 312 > input_v_wavefunction_filename) 313 end if 314 end if 315 316 if(.not.btdb_cget(rtdb,'cpmd:output_v_wavefunction_filename', 317 > 1,output_v_wavefunction_filename)) then 318 if(.not.btdb_cget(rtdb,'pspw:output vvectors', 319 > 1,output_v_wavefunction_filename)) then 320 call util_file_prefix('vmovecs', 321 > output_v_wavefunction_filename) 322 end if 323 end if 324 325 if(.not.btdb_cget(rtdb,'cpmd:xyz_filename', 326 > 1,xyz_filename)) 327 > call util_file_prefix('xyz',xyz_filename) 328! ind = index(exchange_correlation,' ') - 1 329 if(.not.btdb_cget(rtdb,'cpmd:exchange_correlation', 330 > 1,exchange_correlation)) 331 > exchange_correlation = 'vosko' 332 333!$OMP single 334 if (nwpwxc_rtdb_load(rtdb,"dft")) then 335c call nwpwxc_print() 336 has_disp = nwpwxc_has_disp() 337 endif 338!$OMP end single 339 340#include "control_gga.fh" 341 342 343 if(.not. btdb_get(rtdb,'cpmd:geometry_optimize',mt_log,1,move)) 344 > move = .true. 345 if(.not. btdb_get(rtdb,'cpmd:fractional_coordinates', 346 > mt_log,1,frac_coord)) 347 > frac_coord = .false. 348 if(.not. btdb_get(rtdb,'cpmd:npsp',mt_int,1,npsp)) 349 > npsp = 0 350 if(.not. btdb_get(rtdb,'cpmd:fake_mass',mt_dbl,1,fake_mass)) 351 > fake_mass = 800.0d0 352 if(.not. btdb_get(rtdb,'cpmd:time_step',mt_dbl,1,time_step)) 353 > time_step = 5.0d0 354 if(.not. btdb_get(rtdb,'cpmd:loop',mt_int,2,loop)) then 355 loop(1) = 10 356 loop(2) = 100 357 end if 358 if(.not. btdb_get(rtdb,'cpmd:scaling',mt_dbl,2,scaling)) then 359 scaling(1) = 1.0d0 360 scaling(2) = 1.0d0 361 end if 362 363 tolerances(1) = 0.0d0 364 tolerances(2) = 0.0d0 365 tolerances(3) = 0.0d0 366 if(.not. btdb_get(rtdb,'cpmd:ecut',mt_dbl,1,ecut)) 367 > ecut=9000.0d0 368 if(.not. btdb_get(rtdb,'cpmd:wcut',mt_dbl,1,wcut)) 369 > wcut = ecut 370 if(.not. btdb_get(rtdb,'cpmd:rcut',mt_dbl,1,rcut)) 371 > rcut = 0.0d0 372 if(.not. btdb_get(rtdb,'cpmd:ncut',mt_int,1,ncut)) 373 > ncut = 1 374 SA = .true. 375 if (.not.btdb_get(rtdb,'cpmd:sa_decay',mt_dbl,2,sa_decay)) then 376 SA = .false. 377 sa_decay(1) = 1.0d0 378 sa_decay(2) = 1.0d0 379 end if 380 381 if (.not.btdb_get(rtdb,'nwpw:dipole_motion',mt_log, 382 > 1,dipole_motion)) 383 > dipole_motion = .false. 384 385 if (.not.btdb_get(rtdb,'cpmd:fei',mt_log,1,fei)) 386 > fei = .false. 387 388 if (.not.btdb_get(rtdb,'cpmd:fei_quench',mt_log,1,fei_quench)) 389 > fei_quench = .false. 390 391 value3 = .false. 392 if (.not.btdb_get(rtdb,'nwpw:dof_rotation',mt_log,1,dof_rotation)) 393 > then 394 dof_rotation = .false. 395 value3= .true. 396 end if 397 398 if (.not.btdb_get(rtdb,'nwpw:rotation',mt_log,1,rotation)) then 399 rotation = .true. 400 else 401 if (value3) dof_rotation = rotation 402 end if 403 404 405* **** get thermostat information **** 406 if (.not.btdb_get(rtdb,'cpmd:nose',mt_log,1,nose)) 407 > nose = .false. 408 if (.not.btdb_get(rtdb,'cpmd:Pe',mt_dbl,1,Pe)) 409 > Pe = 1200.0d0 410 if (.not.btdb_get(rtdb,'cpmd:Te',mt_dbl,1,Te)) 411 > Te = 298.15d0 412 if (.not.btdb_get(rtdb,'cpmd:Pr',mt_dbl,1,Pr)) 413 > Pr = 1200.0d0 414 if (.not.btdb_get(rtdb,'cpmd:Tr',mt_dbl,1,Tr)) 415 > Tr = 298.15d0 416 417 if (.not.btdb_get(rtdb,'nwpw:spin_orbit',mt_log,1,spin_orbit)) 418 > spin_orbit=.false. 419 if (spin_orbit) ispin=2 420 421 422* ************************************ 423* **** cgsd: stuff or paw: stuff **** 424* ************************************ 425 else if ((code.eq.3).or.(code.eq.8).or. 426 > (code.eq.11).or.(code.eq.12).or.(code.eq.15)) then 427 if (.not.btdb_cget(rtdb,'cgsd:cell_name',1,cell_name)) then 428 cell_name = 'cell_default' 429 end if 430c 431c **** Figure input/output MO vectors **** 432c 433 if (.not. btdb_cget(rtdb, 'pspw:input vectors', 434 > 1,input_wavefunction_filename)) then 435 input_wavefunction_filename = 'atomic' 436 end if 437 438 if (.not. btdb_cget(rtdb, 'pspw:output vectors', 439 > 1,output_wavefunction_filename)) then 440 output_wavefunction_filename = ' ' 441 442 if (output_wavefunction_filename.eq.' ')then 443 if (input_wavefunction_filename.eq.'atomic')then 444 call util_file_prefix('movecs', 445 > output_wavefunction_filename) 446 else 447 output_wavefunction_filename=input_wavefunction_filename 448 endif 449 endif 450 if (input_wavefunction_filename.eq.'atomic')then 451 input_wavefunction_filename = output_wavefunction_filename 452 end if 453 end if 454 455 if (.not.btdb_cget(rtdb,'cgsd:input_ewavefunction_filename', 456 > 1,input_ewavefunction_filename)) 457 > call util_file_prefix('emovecs',input_ewavefunction_filename) 458 459 if (.not.btdb_cget(rtdb,'cgsd:output_ewavefunction_filename', 460 > 1,output_ewavefunction_filename)) 461 > call util_file_prefix('emovecs',output_ewavefunction_filename) 462 463 if (code.eq.11) then 464 465 if(.not.btdb_cget(rtdb,'pspw:input vvectors', 466 > 1,input_v_wavefunction_filename)) then 467 if(.not.btdb_cget(rtdb,'cpmd:input_v_wavefunction_filename', 468 > 1,input_v_wavefunction_filename)) then 469 call util_file_prefix('vmovecs', 470 > input_v_wavefunction_filename) 471 end if 472 end if 473 474 if(.not.btdb_cget(rtdb,'pspw:output vvectors', 475 > 1,output_v_wavefunction_filename)) then 476 if(.not.btdb_cget(rtdb,'cpmd:output_v_wavefunction_filename', 477 > 1,output_v_wavefunction_filename)) then 478 479 call util_file_prefix('vmovecs', 480 > output_v_wavefunction_filename) 481 end if 482 end if 483 end if 484 485 486 if (.not.btdb_cget(rtdb,'cgsd:exchange_correlation', 487 > 1,exchange_correlation)) 488 > exchange_correlation = 'vosko' 489 490!$OMP single 491 if (nwpwxc_rtdb_load(rtdb,"dft")) then 492c call nwpwxc_print() 493 has_disp = nwpwxc_has_disp() 494 endif 495!$OMP end single 496 497#include "control_gga.fh" 498 499* ***** motion filenames **** 500 if(.not.btdb_cget(rtdb,'nwpw:xyz_filename', 501 > 1,xyz_filename)) 502 > call util_file_prefix('xyz',xyz_filename) 503 504 505* **** set Kohn-Sham scf parameters *** 506 if (.not. btdb_get(rtdb,'nwpw:ks_alpha',mt_dbl,1,ks_alpha)) 507 > ks_alpha = 0.25d0 508 if (.not. btdb_get(rtdb,'nwpw:fractional_alpha', 509 > mt_dbl,1,fractional_alpha)) 510 > fractional_alpha = 1.00d0 511 if (.not.btdb_get(rtdb,'nwpw:scf_algorithm', 512 > mt_int,1,scf_algorithm)) 513 > scf_algorithm = 3 514 515 if (.not.btdb_get(rtdb,'nwpw:ks_algorithm', 516 > mt_int,1,ks_algorithm)) 517 > ks_algorithm = 0 518 if (.not.btdb_get(rtdb,'nwpw:kerker_g0', 519 > mt_dbl,1,kerker_g0)) 520 > kerker_g0 = 0.0d0 521 522* **** set maxit_orb maxit_orbs *** 523 if (.not.btdb_get(rtdb, 524 > 'nwpw:ks_maxit_orb',mt_int,1,maxit_orb)) 525 > maxit_orb = 5 526 if (.not.btdb_get(rtdb, 527 > 'nwpw:ks_maxit_orbs',mt_int,1,maxit_orbs)) 528 > maxit_orbs = 1 529 530 531 if (.not.btdb_get(rtdb,'cgsd:npsp',mt_int,1,npsp)) 532 > npsp = 0 533 if (.not.btdb_get(rtdb,'cgsd:fake_mass',mt_dbl,1,fake_mass)) 534 > fake_mass = 400000.0d0 535 if (.not.btdb_get(rtdb,'cgsd:time_step',mt_dbl,1,time_step)) 536 > time_step = 5.8d0 537 if (.not.btdb_get(rtdb,'cgsd:loop',mt_int,2,loop)) then 538 loop(1) = 10 539 loop(2) = 100 540 end if 541 if (.not.btdb_get(rtdb,'cgsd:tolerances',mt_dbl,3,tolerances))then 542 tolerances(1) = 1.0d-7 543 tolerances(2) = 1.0d-7 544 tolerances(3) = 1.0d-4 545 end if 546 547 if(.not. btdb_get(rtdb,'nwpw:scaling',mt_dbl,2,scaling)) then 548 scaling(1) = 1.0d0 549 scaling(2) = 1.0d0 550 end if 551 552 if (.not.btdb_get(rtdb,'cgsd:fractional_coordinates', 553 > mt_log,1,frac_coord)) 554 > frac_coord = .false. 555 556 if (.not.btdb_get(rtdb,'cgsd:ecut',mt_dbl,1,ecut)) 557 > ecut=9000.0d0 558 if (.not.btdb_get(rtdb,'cgsd:wcut',mt_dbl,1,wcut)) 559 > wcut = ecut 560 if (.not.btdb_get(rtdb,'cgsd:rcut',mt_dbl,1,rcut)) 561 > rcut = 0.0d0 562 if (.not.btdb_get(rtdb,'cgsd:ncut',mt_int,1,ncut)) 563 > ncut = 1 564 if (.not.btdb_get(rtdb,'cgsd:mult',mt_int,1,multiplicity)) 565 > multiplicity = 1 566 if (.not.btdb_get(rtdb,'cgsd:ispin',mt_int,1,ispin)) 567 > ispin = 1 568 569 570* **** BO parameterss *** 571 if (.not.btdb_get(rtdb,'nwpw:bo_steps',mt_int,2,bo_steps)) then 572 bo_steps(1) = 10 573 bo_steps(2) = 100 574 end if 575 if (.not.btdb_get(rtdb,'nwpw:bo_time_step',mt_dbl,1,bo_time_step)) 576 > bo_time_step = time_step 577 if (.not.btdb_get(rtdb,'nwpw:bo_algorithm',mt_int,1,bo_algorithm)) 578 > bo_algorithm = 0 579 if (.not.btdb_get(rtdb,'nwpw:bo_fake_mass',mt_dbl,1,bo_fake_mass)) 580 > bo_fake_mass = 500.0d0 581 582 value3 = .false. 583 if (.not.btdb_get(rtdb,'nwpw:dof_rotation',mt_log,1,dof_rotation)) 584 > then 585 dof_rotation = .false. 586 value3= .true. 587 end if 588 589 if (.not.btdb_get(rtdb,'nwpw:rotation',mt_log,1,rotation)) then 590 rotation = .true. 591 else 592 if (value3) dof_rotation = rotation 593 end if 594 595 596 SA = .true. 597 if (.not.btdb_get(rtdb,'nwpw:sa_decay',mt_dbl,2,sa_decay)) then 598 SA = .false. 599 sa_decay(1) = 1.0d0 600 sa_decay(2) = 1.0d0 601 end if 602 603 if (.not.btdb_get(rtdb,'nwpw:dipole_motion',mt_log, 604 > 1,dipole_motion)) 605 > dipole_motion = .false. 606 607 if (.not.btdb_get(rtdb,'nwpw:fei',mt_log,1,fei)) 608 > fei = .false. 609 610 if (.not.btdb_get(rtdb,'nwpw:fei_quench',mt_log,1,fei_quench)) 611 > fei_quench = .false. 612 613 614 615* **** get thermostat information **** 616 if (.not.btdb_get(rtdb,'nwpw:nose',mt_log,1,nose)) 617 > nose = .false. 618 if (.not.btdb_get(rtdb,'nwpw:Pe',mt_dbl,1,Pe)) 619 > Pe = 1200.0d0 620 if (.not.btdb_get(rtdb,'nwpw:Te',mt_dbl,1,Te)) 621 > Te = 298.15d0 622 if (.not.btdb_get(rtdb,'nwpw:Pr',mt_dbl,1,Pr)) 623 > Pr = 1200.0d0 624 if (.not.btdb_get(rtdb,'nwpw:Tr',mt_dbl,1,Tr)) 625 > Tr = 298.15d0 626 627 628* *************************** 629* **** pspw_dplot: stuff **** 630* *************************** 631 else if (code.eq.4) then 632 if (.not.btdb_cget(rtdb,'cgsd:cell_name',1,cell_name)) then 633 cell_name = 'cell_default' 634 end if 635 value = .true. 636 if (.not.btdb_cget(rtdb,'pspw_dplot:wavefunction_filename', 637 > 1,input_wavefunction_filename)) 638 > call util_file_prefix('movecs',input_wavefunction_filename) 639 call psi_get_header(i,ngrid,unita,ispin0,ne) 640 call dcopy(9,unita,1,unita_frozen,1) 641 if (i.eq.3) boundry = 'periodic' 642 if (i.eq.4) boundry = 'aperiodic' 643 644* **** dummy variables **** 645 move = .false. 646 frac_coord = .false. 647 gga = 0 648 fake_mass = 400000.0d0 649 time_step = 5.8d0 650 loop(1) = 0 651 loop(2) = 0 652 tolerances(1) = 1.0d-9 653 tolerances(2) = 1.0d-9 654 tolerances(3) = 1.0d-4 655 if(.not.btdb_get(rtdb,'cgsd:ecut',mt_dbl,1,ecut)) ecut=9000.0d0 656 if(.not.btdb_get(rtdb,'cgsd:wcut',mt_dbl,1,wcut)) wcut=ecut 657 rcut = 0.0d0 658 ncut = 1 659 npsp = 0 660 661c control_read = value 662c return 663 664* ********************* 665* **** band: stuff **** 666* ********************* 667 else if (code.eq.5) then 668 669 if (.not.btdb_cget(rtdb,'band:cell_name',1,cell_name)) then 670 cell_name = 'cell_default' 671 end if 672 673c **** Figure input/output MO vectors **** 674c 675 if (.not.btdb_cget(rtdb,'pspw:input vectors', 676 > 1,input_wavefunction_filename)) then 677 input_wavefunction_filename = 'atomic' 678 end if 679 680 if (.not.btdb_cget(rtdb, 'pspw:output vectors', 681 > 1,output_wavefunction_filename)) then 682 output_wavefunction_filename = ' ' 683 684 if (output_wavefunction_filename.eq.' ')then 685 if (input_wavefunction_filename.eq.'atomic') then 686 call util_file_prefix('movecs', 687 > output_wavefunction_filename) 688 else 689 output_wavefunction_filename=input_wavefunction_filename 690 endif 691 end if 692 end if 693 694 if (input_wavefunction_filename.eq.'atomic')then 695 input_wavefunction_filename = output_wavefunction_filename 696 end if 697 698 if (.not.btdb_cget(rtdb,'cgsd:input_ewavefunction_filename', 699 > 1,input_ewavefunction_filename)) 700 > call util_file_prefix('emovecs',input_ewavefunction_filename) 701 702 if (.not.btdb_cget(rtdb,'cgsd:output_ewavefunction_filename', 703 > 1,output_ewavefunction_filename)) 704 > call util_file_prefix('emovecs',output_ewavefunction_filename) 705 706 707 if (.not.btdb_cget(rtdb,'band:exchange_correlation', 708 > 1,exchange_correlation)) 709 > exchange_correlation = 'vosko' 710 711!$OMP single 712 if (nwpwxc_rtdb_load(rtdb,"dft")) then 713c call nwpwxc_print() 714 has_disp = nwpwxc_has_disp() 715 endif 716!$OMP end single 717 718#include "control_gga.fh" 719 720 721 if(.not.btdb_get(rtdb,'band:geometry_optimize',mt_log,1,move)) 722 > move = .false. 723 if (.not.btdb_get(rtdb,'band:fake_mass',mt_dbl,1,fake_mass)) 724 > fake_mass = 400000.0d0 725 if (.not.btdb_get(rtdb,'band:time_step',mt_dbl,1,time_step)) 726 > time_step = 5.8d0 727 if (.not.btdb_get(rtdb,'band:loop',mt_int,2,loop)) then 728 loop(1) = 10 729 loop(2) = 100 730 end if 731 if(.not.btdb_get(rtdb,'band:tolerances',mt_dbl,3,tolerances)) then 732 tolerances(1) = 1.0d-7 733 tolerances(2) = 1.0d-7 734 tolerances(3) = 1.0d-4 735 end if 736 scaling(1) = 0.0d0 737 scaling(2) = 0.0d0 738 739 if (.not.btdb_get(rtdb,'band:ecut',mt_dbl,1,ecut)) 740 > ecut = 9000.0d0 741 if (.not.btdb_get(rtdb,'band:wcut',mt_dbl,1,wcut)) 742 > wcut = ecut 743 if (.not.btdb_get(rtdb,'band:rcut',mt_dbl,1,rcut)) 744 > rcut = 0.0d0 745 if (.not.btdb_get(rtdb,'band:ncut',mt_int,1,ncut)) 746 > ncut = 1 747 if (.not.btdb_get(rtdb,'band:mult',mt_int,1,multiplicity)) 748 > multiplicity = 1 749 if (.not.btdb_get(rtdb,'band:ispin',mt_int,1,ispin)) 750 > ispin = 1 751 752 if (.not.btdb_get(rtdb,'band:spin_orbit',mt_log,1,spin_orbit)) 753 > spin_orbit=.false. 754 if (spin_orbit) ispin=2 755 756 757 758* **** set Kohn-Sham scf parameters *** 759 if (.not. btdb_get(rtdb,'nwpw:ks_alpha',mt_dbl,1,ks_alpha)) 760 > ks_alpha = 0.25d0 761 if (.not. btdb_get(rtdb,'nwpw:fractional_alpha', 762 > mt_dbl,1,fractional_alpha)) 763 > fractional_alpha = 1.00d0 764 if (.not.btdb_get(rtdb,'nwpw:scf_algorithm', 765 > mt_int,1,scf_algorithm)) 766 > scf_algorithm = 3 767 if (.not.btdb_get(rtdb,'nwpw:ks_algorithm', 768 > mt_int,1,ks_algorithm)) 769 > ks_algorithm = 0 770 if (.not.btdb_get(rtdb,'nwpw:kerker_g0', 771 > mt_dbl,1,kerker_g0)) 772 > kerker_g0 = 0.0d0 773 774 775* **** set maxit_orb maxit_orbs *** 776 if (.not.btdb_get(rtdb, 777 > 'nwpw:ks_maxit_orb',mt_int,1,maxit_orb)) 778 > maxit_orb = 5 779 if (.not.btdb_get(rtdb, 780 > 'nwpw:ks_maxit_orbs',mt_int,1,maxit_orbs)) 781 > maxit_orbs = 1 782 783 784* *********************** 785* **** paw_sd: stuff **** 786* *********************** 787 else if (code.eq.6) then 788 789 if (.not.btdb_get(rtdb,'cgsd:geometry_optimize',mt_log,1,move)) 790 > move = .false. 791 792 if (.not.btdb_cget(rtdb,'cgsd:cell_name',1,cell_name)) then 793 cell_name = 'cell_default' 794 end if 795 if (.not.btdb_cget(rtdb,'cgsd:input_wavefunction_filename', 796 > 1,input_wavefunction_filename)) 797 > call util_file_prefix('movecs',input_wavefunction_filename) 798 if (.not.btdb_cget(rtdb,'cgsd:output_wavefunction_filename', 799 > 1,output_wavefunction_filename)) 800 > call util_file_prefix('movecs',output_wavefunction_filename) 801 if (.not.btdb_cget(rtdb,'cgsd:exchange_correlation', 802 > 1,exchange_correlation)) 803 > exchange_correlation = 'vosko' 804 805 if (inp_compare(.false.,exchange_correlation,'vosko')) then 806 gga = 0 807 else if (inp_compare(.false.,exchange_correlation,'lda')) then 808 gga = 0 809 else if (inp_compare(.false.,exchange_correlation,'svwn5')) then 810 gga = 0 811 else if (inp_compare(.false.,exchange_correlation,'pbe96')) then 812 gga = 10 813 else if (inp_compare(.false.,exchange_correlation,'blyp')) then 814 gga = 11 815 else if (inp_compare(.false.,exchange_correlation,'revpbe')) then 816 gga = 12 817 else 818 gga = 0 819 end if 820 821!$OMP single 822 if (nwpwxc_rtdb_load(rtdb,"dft")) then 823c call nwpwxc_print() 824 has_disp = nwpwxc_has_disp() 825 endif 826!$OMP end single 827 828 if (.not.btdb_get(rtdb,'cgsd:npsp',mt_int,1,npsp)) 829 > npsp = 0 830 if (.not.btdb_get(rtdb,'cgsd:fake_mass',mt_dbl,1,fake_mass)) 831 > fake_mass = 400000.0d0 832 if (.not.btdb_get(rtdb,'cgsd:time_step',mt_dbl,1,time_step)) 833 > time_step = 5.8d0 834 if (.not.btdb_get(rtdb,'cgsd:loop',mt_int,2,loop)) then 835 loop(1) = 10 836 loop(2) = 100 837 end if 838 if(.not.btdb_get(rtdb,'cgsd:tolerances',mt_dbl,3,tolerances)) then 839 tolerances(1) = 1.0d-7 840 tolerances(2) = 1.0d-7 841 tolerances(3) = 1.0d-4 842 end if 843 scaling(1) = 0.0d0 844 scaling(2) = 0.0d0 845 if (.not.btdb_get(rtdb,'cgsd:ecut',mt_dbl,1,ecut)) 846 > ecut = 9000.0d0 847 if (.not.btdb_get(rtdb,'cgsd:wcut',mt_dbl,1,wcut)) 848 > wcut = ecut 849 if (.not.btdb_get(rtdb,'cgsd:rcut',mt_dbl,1,rcut)) 850 > rcut = 0.0d0 851 if (.not.btdb_get(rtdb,'cgsd:ncut',mt_int,1,ncut)) 852 > ncut = 1 853 if (.not.btdb_get(rtdb,'cgsd:mult',mt_int,1,multiplicity)) 854 > multiplicity = 1 855 if (.not.btdb_get(rtdb,'cgsd:ispin',mt_int,1,ispin)) 856 > ispin = 1 857 858 859* ************************* 860* **** paw_cpmd: stuff **** 861* ************************* 862 else if (code.eq.7) then 863 864 if (.not.btdb_cget(rtdb,'cpmd:cell_name',1,cell_name)) then 865 cell_name = 'cell_default' 866 end if 867 868 if (.not.btdb_cget(rtdb,'cpmd:input_wavefunction_filename', 869 > 1,input_wavefunction_filename)) 870 > call util_file_prefix('movecs',input_wavefunction_filename) 871 if (.not.btdb_cget(rtdb,'cpmd:output_wavefunction_filename', 872 > 1,output_wavefunction_filename)) 873 > call util_file_prefix('movecs',output_wavefunction_filename) 874 if (.not.btdb_cget(rtdb,'cpmd:input_v_wavefunction_filename', 875 > 1,input_v_wavefunction_filename)) 876 > call util_file_prefix('vmovecs',input_v_wavefunction_filename) 877 if (.not.btdb_cget(rtdb,'cpmd:output_v_wavefunction_filename', 878 > 1,output_v_wavefunction_filename)) 879 > call util_file_prefix('vmovecs',output_v_wavefunction_filename) 880 881 882* ***** motion filenames **** 883 if(.not.btdb_cget(rtdb,'cpmd:xyz_filename', 884 > 1,xyz_filename)) 885 > call util_file_prefix('xyz',xyz_filename) 886 887 888 if (.not.btdb_cget(rtdb,'cpmd:exchange_correlation', 889 > 1,exchange_correlation)) 890 > exchange_correlation = 'vosko' 891 892 if (inp_compare(.false.,exchange_correlation,'vosko')) then 893 gga = 0 894 else if (inp_compare(.false.,exchange_correlation,'lda')) then 895 gga = 0 896 else if (inp_compare(.false.,exchange_correlation,'svwn5')) then 897 gga = 0 898 else if (inp_compare(.false.,exchange_correlation,'pbe96')) then 899 gga = 10 900 else if (inp_compare(.false.,exchange_correlation,'blyp')) then 901 gga = 11 902 else if (inp_compare(.false.,exchange_correlation,'revpbe')) then 903 gga = 12 904 else 905 gga = 0 906 end if 907 908!$OMP single 909 if (nwpwxc_rtdb_load(rtdb,"dft")) then 910c call nwpwxc_print() 911 has_disp = nwpwxc_has_disp() 912 endif 913!$OMP end single 914 915 916 if (.not.btdb_get(rtdb,'cpmd:geometry_optimize',mt_log,1,move)) 917 > move = .false. 918 if (.not.btdb_get(rtdb,'cpmd:fractional_coordinates', 919 > mt_log,1,frac_coord)) 920 > frac_coord = .false. 921 if (.not.btdb_get(rtdb,'cpmd:npsp',mt_int,1,npsp)) 922 > npsp = 0 923 if (.not.btdb_get(rtdb,'cpmd:fake_mass',mt_dbl,1,fake_mass)) 924 > fake_mass = 800.0d0 925 if (.not.btdb_get(rtdb,'cpmd:time_step',mt_dbl,1,time_step)) 926 > time_step = 5.0d0 927 if (.not.btdb_get(rtdb,'cpmd:loop',mt_int,2,loop)) then 928 loop(1) = 10 929 loop(2) = 100 930 end if 931 if (.not.btdb_get(rtdb,'cpmd:scaling',mt_dbl,2,scaling)) then 932 scaling(1) = 1.0d0 933 scaling(2) = 1.0d0 934 end if 935 tolerances(1) = 0.0d0 936 tolerances(2) = 0.0d0 937 tolerances(3) = 0.0d0 938 if (.not.btdb_get(rtdb,'cpmd:ecut',mt_dbl,1,ecut)) 939 > ecut = 9000.0d0 940 if (.not.btdb_get(rtdb,'cpmd:wcut',mt_dbl,1,wcut)) 941 > wcut = ecut 942 if (.not.btdb_get(rtdb,'cpmd:rcut',mt_dbl,1,rcut)) 943 > rcut = 0.0d0 944 if (.not.btdb_get(rtdb,'cpmd:ncut',mt_int,1,ncut)) 945 > ncut = 1 946 947 SA = .true. 948 if (.not.btdb_get(rtdb,'cpmd:sa_decay',mt_dbl,2,sa_decay)) then 949 SA = .false. 950 sa_decay(1) = 1.0d0 951 sa_decay(2) = 1.0d0 952 end if 953 954 if (.not.btdb_get(rtdb,'nwpw:dipole_motion',mt_log, 955 > 1,dipole_motion)) 956 > dipole_motion = .false. 957 958 if (.not.btdb_get(rtdb,'cpmd:fei',mt_log,1,fei)) 959 > fei = .false. 960 961 if (.not.btdb_get(rtdb,'cpmd:fei_quench',mt_log,1,fei_quench)) 962 > fei_quench = .false. 963 964 value3 = .false. 965 if (.not.btdb_get(rtdb,'nwpw:dof_rotation',mt_log,1,dof_rotation)) 966 > then 967 dof_rotation = .false. 968 value3= .true. 969 end if 970 971 if (.not.btdb_get(rtdb,'nwpw:rotation',mt_log,1,rotation)) then 972 rotation = .true. 973 else 974 if (value3) dof_rotation = rotation 975 end if 976 977* **** get thermostat information **** 978 if (.not.btdb_get(rtdb,'cpmd:nose',mt_log,1,nose)) 979 > nose = .false. 980 if (.not.btdb_get(rtdb,'cpmd:Pe',mt_dbl,1,Pe)) 981 > Pe = 1200.0d0 982 if (.not.btdb_get(rtdb,'cpmd:Te',mt_dbl,1,Te)) 983 > Te = 298.15d0 984 if (.not.btdb_get(rtdb,'cpmd:Pr',mt_dbl,1,Pr)) 985 > Pr = 1200.0d0 986 if (.not.btdb_get(rtdb,'cpmd:Tr',mt_dbl,1,Tr)) 987 > Tr = 298.15d0 988 989 990* ********************** 991* **** pspw_wannier **** 992* ********************** 993 else if (code.eq.9) then 994 995 if (.not.btdb_cget(rtdb,'cgsd:cell_name',1,cell_name)) then 996 cell_name = 'cell_default' 997 end if 998 999c **** Figure input/output MO vectors **** 1000 1001 if (.not.btdb_cget(rtdb,'wannier:input vectors', 1002 > 1,input_wavefunction_filename)) then 1003 if (.not.btdb_cget(rtdb, 'pspw:input vectors', 1004 > 1,input_wavefunction_filename)) then 1005 input_wavefunction_filename = 'atomic' 1006 end if 1007 end if 1008 1009 if (.not.btdb_cget(rtdb,'wannier:output vectors', 1010 > 1,output_wavefunction_filename)) then 1011 if (.not.btdb_cget(rtdb,'pspw:output vectors', 1012 > 1,output_wavefunction_filename)) then 1013 output_wavefunction_filename = ' ' 1014 end if 1015 end if 1016 1017 if (output_wavefunction_filename.eq.' ')then 1018 if (input_wavefunction_filename.eq.'atomic')then 1019 call util_file_prefix('movecs',output_wavefunction_filename) 1020 else 1021 output_wavefunction_filename = input_wavefunction_filename 1022 endif 1023 endif 1024 1025 if (input_wavefunction_filename.eq.'atomic')then 1026 input_wavefunction_filename = output_wavefunction_filename 1027 end if 1028 1029 call psi_get_header(i,ngrid,unita,ispin0,ne) 1030 call dcopy(9,unita,1,unita_frozen,1) 1031 if (i.eq.3) boundry = 'periodic' 1032 if (i.eq.4) boundry = 'aperiodic' 1033 1034* **** dummy variables **** 1035 move = .false. 1036 frac_coord = .false. 1037 gga = 0 1038 fake_mass = 400000.0d0 1039 time_step = 5.8d0 1040 loop(1) = 0 1041 loop(2) = 0 1042 tolerances(1) = 1.0d-9 1043 tolerances(2) = 1.0d-9 1044 tolerances(3) = 1.0d-4 1045 if(.not.btdb_get(rtdb,'cgsd:ecut',mt_dbl,1,ecut)) ecut=9000.0d0 1046 if(.not.btdb_get(rtdb,'cgsd:wcut',mt_dbl,1,wcut)) wcut=ecut 1047 rcut = 0.0d0 1048 ncut = 0 1049 npsp = 0 1050 1051c control_read = value 1052c return 1053 1054* *************************** 1055* **** band_dplot: stuff **** 1056* *************************** 1057 else if (code.eq.10) then 1058 1059 if (.not.btdb_cget(rtdb,'band:cell_name',1,cell_name)) then 1060 cell_name = 'cell_default' 1061 end if 1062 1063 value = .true. 1064 if (.not.btdb_cget(rtdb,'band_dplot:wavefunction_filename', 1065 > 1,input_wavefunction_filename)) 1066 > call util_file_prefix('movecs',input_wavefunction_filename) 1067 call cpsi_get_header(i,ngrid,unita,ispin0,ne,nbrill) 1068 call dcopy(9,unita,1,unita_frozen,1) 1069 if (i.eq.5) boundry = 'periodic' 1070 1071* **** dummy variables **** 1072 move = .false. 1073 frac_coord = .false. 1074 gga = 0 1075 fake_mass = 400000.0d0 1076 time_step = 5.8d0 1077 loop(1) = 0 1078 loop(2) = 0 1079 tolerances(1) = 1.0d-9 1080 tolerances(2) = 1.0d-9 1081 tolerances(3) = 1.0d-4 1082 if(.not.btdb_get(rtdb,'band:ecut',mt_dbl,1,ecut)) ecut=9000.0d0 1083 if(.not.btdb_get(rtdb,'band:wcut',mt_dbl,1,wcut)) wcut=ecut 1084 rcut = 0.0d0 1085 ncut = 0 1086 npsp = 0 1087 1088c control_read = value 1089c return 1090 1091* ******************************** 1092* **** unknown but dont crash **** 1093* ******************************** 1094 else if (code.gt.99) then 1095 call nwpw_timing_end(50) 1096 control_read = value 1097 return 1098 1099* *************************** 1100* **** unknown code type **** 1101* *************************** 1102 else 1103 value = .false. 1104 write(*,*) "control_read: unknown code type:",code 1105 call nwpw_timing_end(50) 1106 control_read = value 1107 return 1108 end if 1109 1110* ***************************** 1111* ***** symmetry variables **** 1112* ***************************** 1113 if (.not.btdb_get(rtdb,'nwpw:symmetry',mt_int,1,symm_number)) 1114 > symm_number = 0 1115 1116* ********************** 1117* ***** cell: stuff **** 1118* ********************** 1119 l = index(cell_name,' ') - 1 1120 rtdb_unita = cell_name(1:l)//':unita' 1121 rtdb_unitaf = cell_name(1:l)//':unita_frozen' 1122 rtdb_ngrid = cell_name(1:l)//':ngrid' 1123 rtdb_boundry = cell_name(1:l)//':boundry' 1124 rtdb_ngrid_small = cell_name(1:l)//':ngrid_small' 1125 1126 1127* **** define unita and boundary **** 1128 if (.not.btdb_get(rtdb,rtdb_unita,mt_dbl,9,unita)) then 1129 call dcopy(9,0.0d0,0,unita,1) 1130 end if 1131 1132 if (.not.btdb_cget(rtdb,rtdb_boundry,1,boundry)) then 1133 boundry = 'periodic' 1134 end if 1135 call check_unita_for_default(rtdb,unita,rtdb_unita,cell_name) 1136 1137 1138* **** define unita_frozen **** 1139 if (.not.btdb_get(rtdb,'nwpw:frozen_lattice',mt_log,1,frozen)) 1140 > frozen = .true. 1141 1142 if (frozen) then 1143 if (.not.btdb_get(rtdb,rtdb_unitaf,mt_dbl,9,unita_frozen)) then 1144 call dcopy(9,unita,1,unita_frozen,1) 1145 value2 = btdb_parallel(.false.) 1146 if (taskid.eq.MASTER) then 1147 value = value.and. 1148 > btdb_put(rtdb,rtdb_unitaf,mt_dbl,9,unita_frozen) 1149 end if 1150 value2 = btdb_parallel(.true.) 1151 else 1152 if (.not.btdb_get(rtdb,'nwpw:frozen_lattice:thresh', 1153 > mt_dbl,1,thresh)) thresh = 0.05d0 1154 error = 0.0d0 1155 do l=1,3 1156 do i=1,3 1157 error = error + ((unita_frozen(i,l)-unita(i,l)))**2 1158 end do 1159 end do 1160 error = dsqrt(error) 1161 if (error.gt.thresh) then 1162 call dcopy(9,unita,1,unita_frozen,1) 1163 value2 = btdb_parallel(.false.) 1164 if (taskid.eq.MASTER) then 1165 value = value.and. 1166 > btdb_put(rtdb,rtdb_unitaf,mt_dbl,9,unita_frozen) 1167 end if 1168 value2 = btdb_parallel(.true.) 1169 end if 1170 end if 1171 else 1172 call dcopy(9,unita,1,unita_frozen,1) 1173 end if 1174 1175 1176* **** define ngrid based on unita and ecut **** 1177 if (.not.btdb_get(rtdb,rtdb_ngrid,mt_int,3,ngrid)) then 1178 if ((ecut.gt.5000.0d0).and.(wcut.gt.5000.0d0)) then 1179 call control_ecut_wcut_default(rtdb,ecut,wcut) 1180 else if ((ecut.gt.5000.0d0).and.(wcut.lt.5000.0d0)) then 1181 ecut = 2.0d0*wcut 1182 end if 1183 1184 !call control_ngrid_default(rtdb,unita,ecut,mapping,ngrid) 1185 call control_ngrid_default(rtdb,unita_frozen,ecut, 1186 > mapping,ngrid) 1187 !ngrid(1) = 32 1188 !ngrid(2) = 32 1189 !ngrid(3) = 32 1190 end if 1191 1192 1193* **** define ngrid_small **** 1194 has_ngrid_small = .false. 1195 ngrid_small(1) = 0 1196 ngrid_small(2) = 0 1197 ngrid_small(3) = 0 1198 if (btdb_get(rtdb,rtdb_ngrid_small,mt_int,3,ngrid_small)) 1199 > has_ngrid_small = .true. 1200 1201* *** set to false if wannier or dplot code *** 1202 if ((code.eq.4) .or.(code.eq.9).or.(code.eq.10)) 1203 > has_ngrid_small = .false. 1204 1205 1206* **** set fractional (smearing) parameters **** 1207 if (.not. 1208 > btdb_get(rtdb,'nwpw:fractional_orbitals',mt_int,2,frac_ne)) then 1209 frac_ne(1) = 0 1210 frac_ne(2) = 0 1211 end if 1212 fractional = (frac_ne(1).gt.0).or.(frac_ne(2).gt.0) 1213 if (.not.btdb_get(rtdb,'nwpw:fractional_temperature', 1214 > mt_dbl,1,frac_temperature)) then 1215 frac_temperature = 0.0d0 1216 end if 1217 if (.not.btdb_get(rtdb,'nwpw:fractional_smeartype', 1218 > mt_int,1,frac_smeartype)) then 1219 frac_smeartype = 0 1220 end if 1221 1222* **** set attenuation parameter **** 1223 if (.not.btdb_get(rtdb,'nwpw:attenuation',mt_dbl,1,attenuation)) 1224 > attenuation = 0.5d0 1225 1226* **** set out of time variables **** 1227 est_step_time = -1 1228 est_finish_time = -1 1229 call current_second(cpu1_time) 1230 1231* **** set gram_schmidt *** 1232 gram_schmidt = .false. 1233 if (.not.btdb_get(rtdb, 1234 > 'nwpw:gram-schmidt',mt_log,1,gram_schmidt)) 1235 > gram_schmidt = .false. 1236 1237* **** set translation *** 1238 translation = control_allow_translation() 1239 dof_translation = translation 1240 1241* **** set two component pseudopotential 1242 two_comp_ppot = .false. 1243 1244 1245* **** ewald ngrid ************ 1246cccc set to some reasonable default ! 1247 if (.not.btdb_get(rtdb,'nwpw:ewald_ngrid',mt_int,3,ewald_grid)) 1248 > then 1249 ewald_grid(1)=ngrid(1) 1250 ewald_grid(2)=ngrid(2) 1251 ewald_grid(3)=ngrid(3) 1252 end if 1253 1254* **** reset mapping if slab and it doesnot fit **** 1255 if (mapping.eq.1) then 1256 if ((ngrid(2).ne.ngrid(3)).or. 1257 > (ngrid(3).lt.np_dimensions(1))) then 1258 mapping = 2 1259 end if 1260 end if 1261 1262 if (np_default) then 1263 if ( (code.eq.1).or. 1264 > (code.eq.2).or. 1265 > (code.eq.3).or. 1266 > (code.eq.9).or. 1267 > (code.eq.11)) then 1268 nmult(1) = 2 1269 nmult(2) = 3 1270 nmult(3) = 5 1271 do i=1,3 1272 do while ((np_dimensions(1).gt.ngrid(1)).and. 1273 > (mod(np_dimensions(1),nmult(i)).eq.0)) 1274 np_dimensions(1) = np_dimensions(1)/nmult(i) 1275 np_dimensions(2) = np_dimensions(2)*nmult(i) 1276 end do 1277 do while ((np_dimensions(1).gt.ngrid(2)).and. 1278 > (mod(np_dimensions(1),nmult(i)).eq.0)) 1279 np_dimensions(1) = np_dimensions(1)/nmult(i) 1280 np_dimensions(2) = np_dimensions(2)*nmult(i) 1281 end do 1282 do while ((np_dimensions(1).gt.ngrid(3)).and. 1283 > (mod(np_dimensions(1),nmult(i)).eq.0)) 1284 np_dimensions(1) = np_dimensions(1)/nmult(i) 1285 np_dimensions(2) = np_dimensions(2)*nmult(i) 1286 end do 1287 end do 1288 end if 1289 end if 1290 pio = pio.and.(np_dimensions(2).gt.1) 1291 1292* **** minimizer **** 1293 if (.not.btdb_get(rtdb,'nwpw:minimizer',mt_int,1,minimizer)) then 1294 minimizer = 1 ! make the default Grassmann cg 1295 end if 1296 1297 call nwpw_timing_end(50) 1298 control_read = value 1299 return 1300 end 1301* *********************************** 1302* * control_ewald_ngrid() 1303* *********************************** 1304 integer function control_ewald_ngrid(i) 1305 implicit none 1306#include "control.fh" 1307 integer i 1308 control_ewald_ngrid=ewald_grid(i) 1309 return 1310 end 1311* *********************************** 1312* * control_ewald_set_ngrid 1313* *********************************** 1314 subroutine ewald_set_ngrid(enx,eny,enz) 1315 implicit none 1316 integer enx,eny,enz 1317 1318#include "control.fh" 1319 1320 ewald_grid(1)=enx 1321 ewald_grid(2)=eny 1322 ewald_grid(3)=enz 1323 return 1324 end 1325 1326* *********************************** 1327* * * 1328* * control_ngrid_default * 1329* * * 1330* *********************************** 1331 subroutine control_ngrid_default(rtdb,unita,ecut,mapping,ngrid) 1332 implicit none 1333 integer rtdb 1334 real*8 unita(3,3),ecut 1335 integer mapping 1336 integer ngrid(3) 1337 1338#include "bafdecls.fh" 1339#include "btdb.fh" 1340#include "errquit.fh" 1341 1342* **** local variables **** 1343 real*8 unitg(3,3),omega 1344 real*8 gx,gy,gz 1345 real*8 xh,yh,zh 1346 integer control_set_ngrid 1347 external control_set_ngrid 1348 1349 call get_cube(unita,unitg,omega) 1350 1351 gx = unitg(1,1) 1352 gy = unitg(2,1) 1353 gz = unitg(3,1) 1354 xh = dsqrt(2.0d0*ecut/(gx*gx + gy*gy + gz*gz))+0.5d0 1355 1356 gx = unitg(1,2) 1357 gy = unitg(2,2) 1358 gz = unitg(3,2) 1359 yh = dsqrt(2.0d0*ecut/(gx*gx + gy*gy + gz*gz))+0.5d0 1360 1361 gx = unitg(1,3) 1362 gy = unitg(2,3) 1363 gz = unitg(3,3) 1364 zh = dsqrt(2.0d0*ecut/(gx*gx + gy*gy + gz*gz))+0.5d0 1365 1366 if (mapping.ge.2) then 1367c ngrid(1) = control_set_ngrid(2.0d0*xh,.true.) 1368c ngrid(2) = control_set_ngrid(2.0d0*yh,.false.) 1369c ngrid(3) = control_set_ngrid(2.0d0*zh,.false.) 1370 ngrid(1) = control_set_ngrid(2.0d0*xh,.true.) 1371 ngrid(2) = control_set_ngrid(2.0d0*yh,.true.) 1372 ngrid(3) = control_set_ngrid(2.0d0*zh,.true.) 1373 else 1374 ngrid(1) = control_set_ngrid(2.0d0*xh,.true.) 1375 ngrid(2) = control_set_ngrid(2.0d0*yh,.true.) 1376 ngrid(3) = control_set_ngrid(2.0d0*zh,.true.) 1377 if (ngrid(2).gt.ngrid(3)) then 1378 ngrid(3) = ngrid(2) 1379 else 1380 ngrid(2) = ngrid(3) 1381 end if 1382 end if 1383 1384 1385c* *** write unita to rtdb - should happen only once during a simulation *** 1386c if (.not.btdb_put(rtdb,rtdb_ngrid,mt_int,3,ngrid)) then 1387c call errquit('cannot write ngrid to rtdb',0,0) 1388c end if 1389 1390 return 1391 end 1392 1393* *********************************** 1394* * * 1395* * control_set_ngrid * 1396* * * 1397* *********************************** 1398* 1399* return n so that it is a multiple of 2,3,5,7 1400* 1401 integer function control_set_ngrid(x,mult2) 1402 implicit none 1403 real*8 x 1404 logical mult2 1405 1406* **** local variables **** 1407 integer nx,ntest 1408 integer nf2,nf3,nf5,nf7 1409 1410 integer factor_count2 1411 external factor_count2 1412 1413* **** find prime factors of nx *** 1414 nx = (x+0.5d0) !*** crude rounding 1415 if ((mult2).and.(mod(nx,2).ne.0)) nx = nx+1 1416 1417 nf2 = factor_count2(nx,2) 1418 nf3 = factor_count2(nx,3) 1419 nf5 = factor_count2(nx,5) 1420 nf7 = factor_count2(nx,7) 1421 ntest = (2**nf2) * (3**nf3) * (5**nf5) * (7**nf7) 1422 do while (nx .ne. ntest) 1423 nx = nx + 1 1424 if (mult2) nx = nx + 1 1425 nf2 = factor_count2(nx,2) 1426 nf3 = factor_count2(nx,3) 1427 nf5 = factor_count2(nx,5) 1428 nf7 = factor_count2(nx,7) 1429 ntest = (2**nf2) * (3**nf3) * (5**nf5) * (7**nf7) 1430 end do 1431 1432 control_set_ngrid = nx 1433 return 1434 end 1435 1436 integer function factor_count2(n,m) 1437 implicit none 1438 integer n,m 1439 integer f,nn 1440 1441 f = 0 1442 nn = n 1443 1444 do while (mod(nn,m).eq.0) 1445 nn = nn/m 1446 f = f + 1 1447 end do 1448 1449 factor_count2 = f 1450 return 1451 end 1452 1453 1454 1455 1456 1457* *********************************** 1458* * * 1459* * check_unita_for_default * 1460* * * 1461* *********************************** 1462 subroutine check_unita_for_default(rtdb,unita,rtdb_unita, 1463 > cell_name) 1464 implicit none 1465 integer rtdb 1466 real*8 unita(3,3) 1467 character*50 rtdb_unita 1468 character*50 cell_name 1469 1470 1471 1472#include "bafdecls.fh" 1473#include "btdb.fh" 1474#include "beom.fh" 1475#include "errquit.fh" 1476 1477* **** local variables **** 1478 integer MASTER,taskid 1479 parameter (MASTER=0) 1480 logical value,box_orient 1481 integer geom,box_type,l,isystype 1482 real*8 box_delta 1483 character*50 rtdb_name 1484 1485 value = (unita(1,1) .eq. 0.0d0).and. 1486 > (unita(2,1) .eq. 0.0d0).and. 1487 > (unita(3,1) .eq. 0.0d0).and. 1488 > (unita(1,2) .eq. 0.0d0).and. 1489 > (unita(2,2) .eq. 0.0d0).and. 1490 > (unita(3,2) .eq. 0.0d0).and. 1491 > (unita(1,3) .eq. 0.0d0).and. 1492 > (unita(2,3) .eq. 0.0d0).and. 1493 > (unita(3,3) .eq. 0.0d0) 1494 1495 if (value) then 1496 value = beom_create(geom,'geometry') 1497 value = value.and.beom_rtdb_load(rtdb,geom,'geometry') 1498 value = value.and.geom_amatrix_get(geom,unita) 1499 value = value.and.geom_systype_get(geom,isystype) 1500 if (.not. value) call errquit('cannot load geometry',0, 1501 & GEOM_ERR) 1502 1503 value = (unita(1,1) .eq. 1.0d0).and. 1504 > (unita(2,1) .eq. 0.0d0).and. 1505 > (unita(3,1) .eq. 0.0d0).and. 1506 > (unita(1,2) .eq. 0.0d0).and. 1507 > (unita(2,2) .eq. 1.0d0).and. 1508 > (unita(3,2) .eq. 0.0d0).and. 1509 > (unita(1,3) .eq. 0.0d0).and. 1510 > (unita(2,3) .eq. 0.0d0).and. 1511 > (unita(3,3) .eq. 1.0d0) 1512 if (value) then 1513 l = index(cell_name,' ') - 1 1514 rtdb_name = cell_name(1:l)//':box_delta' 1515 if (.not.btdb_get(rtdb,rtdb_name,mt_dbl,1,box_delta)) 1516 > box_delta = 5.0d0 1517 rtdb_name = cell_name(1:l)//':box_type' 1518 if (.not.btdb_get(rtdb,rtdb_name,mt_int,1,box_type)) 1519 > box_type = 1 1520 rtdb_name = cell_name(1:l)//':box_orient' 1521 if (.not.btdb_get(rtdb,rtdb_name,mt_log,1, 1522 > box_orient)) 1523 > box_orient = .true. 1524 call control_find_box(geom,box_type,box_orient,box_delta, 1525 > unita) 1526 1527* *** write unita to rtdb - should happen only once during a simulation *** 1528 call Parallel_taskid(taskid) 1529 value = btdb_parallel(.false.) 1530 if (taskid.eq.MASTER) then 1531 if (.not.btdb_put(rtdb,rtdb_unita,mt_dbl,9,unita)) then 1532 call errquit('cannot write unita to rtdb',0,0) 1533 end if 1534 end if 1535 value = btdb_parallel(.true.) 1536 1537 !unita(1,1) = 20.0d0 1538 !unita(2,1) = 0.0d0 1539 !unita(3,1) = 0.0d0 1540 !unita(1,2) = 0.0d0 1541 !unita(2,2) = 20.0d0 1542 !unita(3,2) = 0.0d0 1543 !unita(1,3) = 0.0d0 1544 !unita(2,3) = 0.0d0 1545 !unita(3,3) = 20.0d0 1546 end if 1547 value = beom_destroy(geom) 1548 if (.not. value) call errquit('cannot destroy geom',0, 1549 & GEOM_ERR) 1550 1551 end if 1552 1553 return 1554 end 1555 1556* *********************************** 1557* * * 1558* * control_find_box * 1559* * * 1560* *********************************** 1561 subroutine control_find_box(geom,box_type,orient,delta,unita) 1562 implicit none 1563 integer geom 1564 integer box_type 1565 logical orient 1566 real*8 delta 1567 real*8 unita(3,3) 1568 1569#include "bafdecls.fh" 1570#include "btdb.fh" 1571#include "beom.fh" 1572#include "errquit.fh" 1573 1574* *** local variables *** 1575 integer lwork 1576 parameter (lwork=9) 1577 double precision work(lwork) 1578 1579 integer ii,nion,ierr,i,j 1580 double precision q,rxyz(3),mtensor(3,3),frac(3),L1,L2,L3 1581 double precision x,y,z,meig(3) 1582 double precision Ixx,Iyx,Izx 1583 double precision Ixy,Iyy,Izy 1584 double precision Ixz,Iyz,Izz,Lmax,tmass 1585 double precision txy,txz,tyz,txx_yy,txx_zz,tyy_zz,tx,ty,tz 1586 1587 character*16 t 1588 logical value 1589 1590* **** external functions **** 1591 logical control_notqmmmq 1592 external control_notqmmmq 1593 1594 if (.not.geom_ncent(geom,nion)) then 1595 call errquit('cannot load nion from geom',0,GEOM_ERR) 1596 end if 1597 1598* *** find principle axes wrt origin (0,0,0) *** 1599 mtensor(1,1) = 1.0d0 1600 mtensor(2,1) = 0.0d0 1601 mtensor(3,1) = 0.0d0 1602 mtensor(1,2) = 0.0d0 1603 mtensor(2,2) = 1.0d0 1604 mtensor(3,2) = 0.0d0 1605 mtensor(1,3) = 0.0d0 1606 mtensor(2,3) = 0.0d0 1607 mtensor(3,3) = 1.0d0 1608 if (orient) then 1609 call dcopy(9,0.0d0,0,mtensor,1) 1610 txy = 0.0d0 1611 txz = 0.0d0 1612 tyz = 0.0d0 1613 txx_yy = 0.0d0 1614 txx_zz = 0.0d0 1615 tyy_zz = 0.0d0 1616 tx = 0.0d0 1617 ty = 0.0d0 1618 tz = 0.0d0 1619 do ii=1,nion 1620 value = geom_cent_get(geom,ii,t,rxyz,q) 1621 if (.not.control_notqmmmq(t)) then 1622 x = rxyz(1) 1623 y = rxyz(2) 1624 z = rxyz(3) 1625 1626 txy = txy + x*y 1627 txz = txz + x*z 1628 tyz = tyz + y*z 1629 1630 txx_yy = txx_yy + (x*x + y*y) 1631 txx_zz = txx_zz + (x*x + z*z) 1632 tyy_zz = tyy_zz + (y*y + z*z) 1633 1634 tx = tx + x 1635 ty = ty + y 1636 tz = tz + z 1637 end if 1638 end do 1639 tmass = dble(nion) 1640 Ixx = tyy_zz - ty*ty/tmass - tz*tz/tmass 1641 Iyy = txx_zz - tx*tx/tmass - tz*tz/tmass 1642 Izz = txx_yy - tx*tx/tmass - ty*ty/tmass 1643 Ixy = txy - tx*ty/tmass 1644 Iyx = Ixy 1645 Ixz = txz - tx*tz/tmass 1646 Izx = Ixz 1647 Iyz = tyz - ty*tz/tmass 1648 Izy = Iyz 1649 1650 mtensor(1,1) = Ixx 1651 mtensor(2,1) = -Ixy 1652 mtensor(3,1) = -Ixz 1653 1654 mtensor(1,2) = -Iyx 1655 mtensor(2,2) = Iyy 1656 mtensor(3,2) = -Iyz 1657 1658 mtensor(1,3) = -Izx 1659 mtensor(2,3) = -Izy 1660 mtensor(3,3) = Izz 1661 1662c **** longest dimension is along a1 **** 1663 ierr = 0 1664 call DSYEV('V','U',3,mtensor,3,meig,work,lwork,ierr) 1665 1666 1667 1668c !*** reorder eigenvectors - make longest dimesion along a3 **** 1669c x = mtensor(1,1) 1670c y = mtensor(2,1) 1671c z = mtensor(3,1) 1672c mtensor(1,1) = mtensor(1,3) 1673c mtensor(2,1) = mtensor(2,3) 1674c mtensor(3,1) = mtensor(3,3) 1675c mtensor(1,3) = x 1676c mtensor(2,3) = y 1677c mtensor(3,3) = z 1678 end if 1679 1680* ***************** 1681* *** cubic box *** 1682* ***************** 1683 if (box_type.eq.0) then 1684 1685* **** define L1 *** 1686 L1 = 0.0d0 1687 do ii=1,nion 1688 if (.not.geom_cent_get(geom,ii,t,rxyz,q)) then 1689 call errquit('cannot load center from geom',0,GEOM_ERR) 1690 end if 1691 if (.not.control_notqmmmq(t)) then 1692 frac(1) = mtensor(1,1)*rxyz(1) 1693 > + mtensor(2,1)*rxyz(2) 1694 > + mtensor(3,1)*rxyz(3) 1695 frac(2) = mtensor(1,2)*rxyz(1) 1696 > + mtensor(2,2)*rxyz(2) 1697 > + mtensor(3,2)*rxyz(3) 1698 frac(3) = mtensor(1,3)*rxyz(1) 1699 > + mtensor(2,3)*rxyz(2) 1700 > + mtensor(3,3)*rxyz(3) 1701 if ((frac(1)+delta).gt.( 0.5*L1)) L1 = 2.0d0*(frac(1)+delta) 1702 if ((frac(2)+delta).gt.( 0.5*L1)) L1 = 2.0d0*(frac(2)+delta) 1703 if ((frac(3)+delta).gt.( 0.5*L1)) L1 = 2.0d0*(frac(3)+delta) 1704 if ((frac(1)-delta).lt.(-0.5*L1)) 1705 > L1 = 2.0d0*(dabs(frac(1))+delta) 1706 if ((frac(2)-delta).lt.(-0.5*L1)) 1707 > L1 = 2.0d0*(dabs(frac(2))+delta) 1708 if ((frac(3)-delta).lt.(-0.5*L1)) 1709 > L1 = 2.0d0*(dabs(frac(3))+delta) 1710 end if 1711 end do 1712 !*** put threshold on smallest cubic box *** 1713 if (L1.lt.26.0d0) L1=26.0d0 1714 1715* **** define unit cell **** 1716 unita(1,1) = L1*mtensor(1,1) 1717 unita(2,1) = L1*mtensor(2,1) 1718 unita(3,1) = L1*mtensor(3,1) 1719 unita(1,2) = L1*mtensor(1,2) 1720 unita(2,2) = L1*mtensor(2,2) 1721 unita(3,2) = L1*mtensor(3,2) 1722 unita(1,3) = L1*mtensor(1,3) 1723 unita(2,3) = L1*mtensor(2,3) 1724 unita(3,3) = L1*mtensor(3,3) 1725 1726* ************************ 1727* *** orthorhombic box *** 1728* ************************ 1729 else if (box_type.eq.1) then 1730 1731* **** define L1, L2, and L3 **** 1732 L1 = 0.0d0 1733 L2 = 0.0d0 1734 L3 = 0.0d0 1735 do ii=1,nion 1736 if (.not.geom_cent_get(geom,ii,t,rxyz,q)) then 1737 call errquit('cannot load center from geom',0,GEOM_ERR) 1738 end if 1739 if (.not.control_notqmmmq(t)) then 1740 frac(1) = mtensor(1,1)*rxyz(1) 1741 > + mtensor(2,1)*rxyz(2) 1742 > + mtensor(3,1)*rxyz(3) 1743 frac(2) = mtensor(1,2)*rxyz(1) 1744 > + mtensor(2,2)*rxyz(2) 1745 > + mtensor(3,2)*rxyz(3) 1746 frac(3) = mtensor(1,3)*rxyz(1) 1747 > + mtensor(2,3)*rxyz(2) 1748 > + mtensor(3,3)*rxyz(3) 1749c write(*,*) "ii,frac=",ii,frac 1750c write(*,*) " =", 1751c > mtensor(1,1),mtensor(2,1),mtensor(3,1) 1752c write(*,*) " =",rxyz 1753c write(*,*) 1754 if ((frac(1)+delta).gt.(0.5d0*L1)) L1 = 2.0d0*(frac(1)+delta) 1755 if ((frac(2)+delta).gt.(0.5d0*L2)) L2 = 2.0d0*(frac(2)+delta) 1756 if ((frac(3)+delta).gt.(0.5d0*L3)) L3 = 2.0d0*(frac(3)+delta) 1757 if ((frac(1)-delta).lt.(-0.5d0*L1)) 1758 > L1 = 2.0d0*(dabs(frac(1))+delta) 1759 if ((frac(2)-delta).lt.(-0.5d0*L2)) 1760 > L2 = 2.0d0*(dabs(frac(2))+delta) 1761 if ((frac(3)-delta).lt.(-0.5d0*L3)) 1762 > L3 = 2.0d0*(dabs(frac(3))+delta) 1763 end if 1764 end do 1765 1766c write(*,*) "L1=",L1 1767c write(*,*) "L2=",L2 1768c write(*,*) "L3=",L3 1769 !*** put threshold on smallest box *** 1770 Lmax = L1 1771 if (Lmax.lt.L2) Lmax = L2 1772 if (Lmax.lt.L3) Lmax = L3 1773 if (Lmax.lt.26.0d0) then 1774 L1 = 19.0d0 1775 L2 = 19.0d0 1776 L3 = 19.0d0 1777 mtensor(1,1) = 1.0d0 1778 mtensor(2,1) = 1.0d0 1779 mtensor(3,1) = 0.0d0 1780 mtensor(1,2) = 1.0d0 1781 mtensor(2,2) = 0.0d0 1782 mtensor(3,2) = 1.0d0 1783 mtensor(1,3) = 0.0d0 1784 mtensor(2,3) = 1.0d0 1785 mtensor(3,3) = 1.0d0 1786 end if 1787 1788* **** define unit cell **** 1789 unita(1,1) = L1*mtensor(1,1) 1790 unita(2,1) = L1*mtensor(2,1) 1791 unita(3,1) = L1*mtensor(3,1) 1792 unita(1,2) = L2*mtensor(1,2) 1793 unita(2,2) = L2*mtensor(2,2) 1794 unita(3,2) = L2*mtensor(3,2) 1795 unita(1,3) = L3*mtensor(1,3) 1796 unita(2,3) = L3*mtensor(2,3) 1797 unita(3,3) = L3*mtensor(3,3) 1798 1799* **** unknown box type **** 1800 else 1801 call errquit('invalid box_type',0,0) 1802 end if 1803 1804 1805 return 1806 end 1807 1808* *************************** 1809* * * 1810* * control_notqmmmq * 1811* * * 1812* *************************** 1813 logical function control_notqmmmq(string) 1814 implicit none 1815 character*16 string 1816 1817 logical qmmmq 1818 1819 qmmmq = .false. 1820 if (index(string,'^').gt.0) qmmmq = .true. 1821 if (index(string,'x').eq.1) qmmmq = .true. 1822 if (index(string,'X').eq.1) qmmmq = .true. 1823 if (index(string,'bq').eq.1) qmmmq = .true. 1824 if (index(string,'Bq').eq.1) qmmmq = .true. 1825 if (index(string,'bQ').eq.1) qmmmq = .true. 1826 if (index(string,'BQ').eq.1) qmmmq = .true. 1827 1828 1829 control_notqmmmq = qmmmq 1830 return 1831 end 1832 1833 1834 1835* *********************************** 1836* * * 1837* * control_move * 1838* * * 1839* *********************************** 1840 logical function control_move() 1841 implicit none 1842 1843#include "control.fh" 1844 1845 control_move = move 1846 return 1847 end 1848 1849 1850* *********************************** 1851* * * 1852* * control_rotation * 1853* * * 1854* *********************************** 1855 logical function control_rotation() 1856 implicit none 1857 1858#include "control.fh" 1859 1860 control_rotation = rotation 1861 return 1862 end 1863 1864* *********************************** 1865* * * 1866* * control_dof_rotation * 1867* * * 1868* *********************************** 1869 logical function control_dof_rotation() 1870 implicit none 1871 1872#include "control.fh" 1873 1874 control_dof_rotation = dof_rotation 1875 return 1876 end 1877 1878* *********************************** 1879* * * 1880* * control_translation * 1881* * * 1882* *********************************** 1883 logical function control_translation() 1884 implicit none 1885 1886#include "control.fh" 1887 1888 control_translation = translation 1889 return 1890 end 1891 1892* *********************************** 1893* * * 1894* * control_dof_translation * 1895* * * 1896* *********************************** 1897 logical function control_dof_translation() 1898 implicit none 1899 1900#include "control.fh" 1901 1902 control_dof_translation = dof_translation 1903 return 1904 end 1905 1906 1907 1908 1909 1910* *********************************** 1911* * * 1912* * control_translate_vector * 1913* * * 1914* *********************************** 1915 subroutine control_translate_vector(rtrans) 1916 implicit none 1917 real*8 rtrans(3) 1918 1919#include "bafdecls.fh" 1920#include "btdb.fh" 1921 1922* **** control_rtdb common block **** 1923 integer rtdb 1924 common / control_rtdb1 / rtdb 1925 1926 if (.not.btdb_get(rtdb,'nwpw:translate_vector', 1927 > mt_dbl,3,rtrans)) then 1928 rtrans(1) = 0.0d0 1929 rtrans(2) = 0.0d0 1930 rtrans(3) = 0.0d0 1931 end if 1932 return 1933 end 1934 1935* *********************************** 1936* * * 1937* * control_translate_reorder * 1938* * * 1939* *********************************** 1940 subroutine control_translate_reorder(reorder) 1941 implicit none 1942 logical reorder 1943 1944#include "bafdecls.fh" 1945#include "btdb.fh" 1946 1947* **** control_rtdb common block **** 1948 integer rtdb 1949 common / control_rtdb1 / rtdb 1950 1951 if (.not.btdb_get(rtdb,'nwpw:translate_reorder', 1952 > mt_log,1,reorder)) then 1953 reorder = .true. 1954 end if 1955 return 1956 end 1957 1958 1959 1960* *********************************** 1961* * * 1962* * control_translate_geom_name * 1963* * * 1964* *********************************** 1965 subroutine control_translate_geom_name(geom_name) 1966 implicit none 1967 character*(*) geom_name 1968 1969#include "bafdecls.fh" 1970#include "btdb.fh" 1971 1972* **** control_rtdb common block **** 1973 integer rtdb 1974 common / control_rtdb1 / rtdb 1975 1976 1977 if (.not.btdb_cget(rtdb,'nwpw:translate_geom_name', 1978 > 1,geom_name)) then 1979 geom_name = "translated_geometry" 1980 end if 1981 return 1982 end 1983 1984 1985 1986 1987* *********************************** 1988* * * 1989* * control_out_of_time * 1990* * * 1991* *********************************** 1992 1993* This function is used to estimate if there is 1994* enough time to perform another iteration. The 1995* routine control_read intializes this routine. To 1996* determine if there is enough time left to do another 1997* iteration this routine uses estimates for the amount 1998* of time to finish the simulation (est_finish_time) and 1999* the amount of time to perform another iteration step 2000* (est_step_time). Where 2001* 2002* est_finish_time = 2*(time elapsed from call to control_read 2003* to the first call to control_out_of_time) 2004* 2005* est_step_time = (time elapsed between successive calls to 2006* control_out_of_time) 2007* 2008* Uses: control_blktime common block located in control.fh 2009* util_test_time_remaining 2010* created: 5-8-2002 2011 2012 logical function control_out_of_time() 2013 implicit none 2014 2015#include "control.fh" 2016 2017* **** control_rtdb common block **** 2018 integer rtdb 2019 common / control_rtdb1 / rtdb 2020 2021* **** local variables **** 2022 logical value 2023 integer required_time 2024 2025* **** external functions **** 2026 logical util_test_time_remaining 2027 external util_test_time_remaining 2028 2029 call current_second(cpu2_time) 2030 2031* **** This is the first time this routine has been called **** 2032 if (est_finish_time.eq.-1) then 2033 est_finish_time = 2*int(cpu2_time-cpu1_time) ! crude estimate 2034 value = .false. 2035 2036* **** This routine has been called two or more times **** 2037 else 2038 est_step_time = int(cpu2_time-cpu1_time)+1 ! no statistical info used 2039 required_time = est_step_time + est_finish_time 2040 value = .not.util_test_time_remaining(rtdb,required_time) 2041 end if 2042 2043 cpu1_time = cpu2_time 2044 2045 control_out_of_time = value 2046 return 2047 end 2048 2049 2050* *********************************** 2051* * * 2052* * control_frac_coord * 2053* * * 2054* *********************************** 2055 logical function control_frac_coord() 2056 implicit none 2057 2058#include "control.fh" 2059 2060 control_frac_coord = frac_coord 2061 return 2062 end 2063 2064 2065 2066 2067 2068* *********************************** 2069* * * 2070* * control_code * 2071* * * 2072* *********************************** 2073 integer function control_code() 2074 implicit none 2075 2076#include "control.fh" 2077 2078 control_code = code 2079 return 2080 end 2081 2082 2083 2084* *********************************** 2085* * * 2086* * control_ngrid * 2087* * * 2088* *********************************** 2089 integer function control_ngrid(ijk) 2090 implicit none 2091 integer ijk 2092 2093#include "control.fh" 2094 2095 control_ngrid = ngrid(ijk) 2096 return 2097 end 2098 2099 2100* *********************************** 2101* * * 2102* * control_ngrid_small * 2103* * * 2104* *********************************** 2105 integer function control_ngrid_small(ijk) 2106 implicit none 2107 integer ijk 2108 2109#include "control.fh" 2110 2111 control_ngrid_small = ngrid_small(ijk) 2112 return 2113 end 2114 2115* *********************************** 2116* * * 2117* * control_has_ngrid_small * 2118* * * 2119* *********************************** 2120 logical function control_has_ngrid_small() 2121 implicit none 2122 2123#include "control.fh" 2124 2125 control_has_ngrid_small = has_ngrid_small 2126 return 2127 end 2128 2129 2130* *********************************** 2131* * * 2132* * control_it_in * 2133* * * 2134* *********************************** 2135 integer function control_it_in() 2136 implicit none 2137 2138#include "control.fh" 2139 2140 control_it_in = loop(1) 2141 return 2142 end 2143 2144 2145* *********************************** 2146* * * 2147* * control_it_out * 2148* * * 2149* *********************************** 2150 integer function control_it_out() 2151 implicit none 2152 2153#include "control.fh" 2154 2155 control_it_out = loop(2) 2156 return 2157 end 2158 2159* *********************************** 2160* * * 2161* * control_bo_steps_in * 2162* * * 2163* *********************************** 2164 integer function control_bo_steps_in() 2165 implicit none 2166 2167#include "control.fh" 2168 2169 control_bo_steps_in = bo_steps(1) 2170 return 2171 end 2172 2173* *********************************** 2174* * * 2175* * control_bo_steps_out * 2176* * * 2177* *********************************** 2178 integer function control_bo_steps_out() 2179 implicit none 2180 2181#include "control.fh" 2182 2183 control_bo_steps_out = bo_steps(2) 2184 return 2185 end 2186 2187* *********************************** 2188* * * 2189* * control_bo_algorithm * 2190* * * 2191* *********************************** 2192 integer function control_bo_algorithm() 2193 implicit none 2194 2195#include "control.fh" 2196 2197 control_bo_algorithm = bo_algorithm 2198 return 2199 end 2200 2201 2202 2203* *********************************** 2204* * * 2205* * control_time_step * 2206* * * 2207* *********************************** 2208 real*8 function control_time_step() 2209 implicit none 2210 2211#include "control.fh" 2212 2213 control_time_step = time_step 2214 return 2215 end 2216 2217* *********************************** 2218* * * 2219* * control_bo_time_step * 2220* * * 2221* *********************************** 2222 real*8 function control_bo_time_step() 2223 implicit none 2224 2225#include "control.fh" 2226 2227 control_bo_time_step = bo_time_step 2228 2229 return 2230 end 2231 2232* *********************************** 2233* * * 2234* * control_ion_time_step * 2235* * * 2236* *********************************** 2237* Used by ion.F, the ion_time_step is 2238* set to time_step if Car-Parrinello 2239* and it is set to bo_time_step if 2240* Born-Oppenheimer. 2241 real*8 function control_ion_time_step() 2242 implicit none 2243 2244#include "control.fh" 2245 2246 !*** BO dynamics **** 2247 if (code.eq.11) then 2248 control_ion_time_step = bo_time_step 2249 2250 !*** CP dynamics **** 2251 else 2252 control_ion_time_step = time_step 2253 end if 2254 2255 return 2256 end 2257 2258 2259 2260* *********************************** 2261* * * 2262* * control_fake_mass * 2263* * * 2264* *********************************** 2265 real*8 function control_fake_mass() 2266 implicit none 2267 2268#include "control.fh" 2269 2270 control_fake_mass = fake_mass 2271 return 2272 end 2273 2274 2275* *********************************** 2276* * * 2277* * control_bo_fake_mass * 2278* * * 2279* *********************************** 2280 real*8 function control_bo_fake_mass() 2281 implicit none 2282 2283#include "control.fh" 2284 2285 control_bo_fake_mass = bo_fake_mass 2286 return 2287 end 2288 2289* *********************************** 2290* * * 2291* * control_ks_alpha * 2292* * * 2293* *********************************** 2294 real*8 function control_ks_alpha() 2295 implicit none 2296 2297#include "control.fh" 2298 2299 control_ks_alpha = ks_alpha 2300 return 2301 end 2302 2303* *********************************** 2304* * * 2305* * control_fractional_alpha * 2306* * * 2307* *********************************** 2308 real*8 function control_fractional_alpha() 2309 implicit none 2310 2311#include "control.fh" 2312 2313 control_fractional_alpha = fractional_alpha 2314 return 2315 end 2316 2317 2318* *********************************** 2319* * * 2320* * control_kerker_g0 * 2321* * * 2322* *********************************** 2323 real*8 function control_kerker_g0() 2324 implicit none 2325 2326#include "control.fh" 2327 2328 control_kerker_g0 = kerker_g0 2329 return 2330 end 2331 2332 2333* *********************************** 2334* * * 2335* * control_scf_algorithm * 2336* * * 2337* *********************************** 2338 integer function control_scf_algorithm() 2339 implicit none 2340 2341#include "control.fh" 2342 2343 control_scf_algorithm = scf_algorithm 2344 return 2345 end 2346 2347* *********************************** 2348* * * 2349* * control_ks_algorithm * 2350* * * 2351* *********************************** 2352 integer function control_ks_algorithm() 2353 implicit none 2354 2355#include "control.fh" 2356 2357 control_ks_algorithm = ks_algorithm 2358 return 2359 end 2360 2361 2362* *********************************** 2363* * * 2364* * control_ks_maxit_orb * 2365* * * 2366* *********************************** 2367 integer function control_ks_maxit_orb() 2368 implicit none 2369 2370#include "control.fh" 2371 2372 control_ks_maxit_orb = maxit_orb 2373 return 2374 end 2375 2376* *********************************** 2377* * * 2378* * control_ks_maxit_orbs * 2379* * * 2380* *********************************** 2381 integer function control_ks_maxit_orbs() 2382 implicit none 2383 2384#include "control.fh" 2385 2386 control_ks_maxit_orbs = maxit_orbs 2387 return 2388 end 2389 2390 2391 2392* *********************************** 2393* * * 2394* * control_tole * 2395* * * 2396* *********************************** 2397 real*8 function control_tole() 2398 implicit none 2399 2400#include "control.fh" 2401 2402 control_tole = tolerances(1) 2403 return 2404 end 2405 2406 2407* *********************************** 2408* * * 2409* * control_tolc * 2410* * * 2411* *********************************** 2412 real*8 function control_tolc() 2413 implicit none 2414 2415#include "control.fh" 2416 2417 control_tolc = tolerances(2) 2418 return 2419 end 2420 2421 2422* *********************************** 2423* * * 2424* * control_tolr * 2425* * * 2426* *********************************** 2427 real*8 function control_tolr() 2428 implicit none 2429 2430#include "control.fh" 2431 2432 control_tolr = tolerances(3) 2433 return 2434 end 2435 2436* *********************************** 2437* * * 2438* * control_rte * 2439* * * 2440* *********************************** 2441 real*8 function control_rte() 2442 implicit none 2443 2444#include "control.fh" 2445 2446 control_rte = scaling(1) 2447 return 2448 end 2449 2450* *********************************** 2451* * * 2452* * control_rti * 2453* * * 2454* *********************************** 2455 real*8 function control_rti() 2456 implicit none 2457 2458#include "control.fh" 2459 2460 control_rti = scaling(2) 2461 return 2462 end 2463 2464 2465* *********************************** 2466* * * 2467* * control_unita * 2468* * * 2469* *********************************** 2470 real*8 function control_unita(i,j) 2471 implicit none 2472 integer i,j 2473 2474#include "control.fh" 2475 2476 control_unita = unita(i,j) 2477 return 2478 end 2479 2480 2481* *********************************** 2482* * * 2483* * control_set_unita * 2484* * * 2485* *********************************** 2486 subroutine control_set_unita(unita_in) 2487 implicit none 2488 real*8 unita_in(*) 2489 2490#include "bafdecls.fh" 2491#include "btdb.fh" 2492#include "errquit.fh" 2493#include "control.fh" 2494 2495* **** local variables **** 2496 integer taskid,MASTER 2497 parameter (MASTER=0) 2498 logical value 2499 integer l 2500 character*50 rtdb_unita 2501 2502* **** control_rtdb common block **** 2503 integer rtdb 2504 common / control_rtdb1 / rtdb 2505 2506 2507 call dcopy(9,unita_in,1,unita,1) 2508 2509* **** put unita on the rtdb **** 2510 l = index(cell_name,' ') - 1 2511 rtdb_unita = cell_name(1:l)//':unita' 2512 2513 call Parallel_taskid(taskid) 2514 value=btdb_parallel(.false.) 2515 if (taskid.eq.MASTER) then 2516 if (.not.btdb_put(rtdb,rtdb_unita,mt_dbl,9,unita)) 2517 > call errquit('control_set_unita:writing unita',0,RTDB_ERR) 2518 end if 2519 value=btdb_parallel(.true.) 2520 return 2521 end 2522 2523 2524* *********************************** 2525* * * 2526* * control_get_unita * 2527* * * 2528* *********************************** 2529 subroutine control_get_unita(unita_out) 2530 implicit none 2531 real*8 unita_out(*) 2532 2533#include "control.fh" 2534 2535 call dcopy(9,unita,1,unita_out,1) 2536 return 2537 end 2538 2539 2540 2541* *********************************** 2542* * * 2543* * control_set_unita_frozen * 2544* * * 2545* *********************************** 2546 subroutine control_set_unita_frozen(unita_in) 2547 implicit none 2548 real*8 unita_in(*) 2549 2550#include "bafdecls.fh" 2551#include "btdb.fh" 2552#include "errquit.fh" 2553#include "control.fh" 2554 2555* **** local variables **** 2556 integer taskid,MASTER 2557 parameter (MASTER=0) 2558 logical value 2559 integer l 2560 character*50 rtdb_unita 2561 2562* **** control_rtdb common block **** 2563 integer rtdb 2564 common / control_rtdb1 / rtdb 2565 2566 call dcopy(9,unita_in,1,unita_frozen,1) 2567 2568* **** put unita_frozen on the rtdb **** 2569 l = index(cell_name,' ') - 1 2570 rtdb_unita = cell_name(1:l)//':unita_frozen' 2571 2572 call Parallel_taskid(taskid) 2573 value=btdb_parallel(.false.) 2574 if (taskid.eq.MASTER) then 2575 if (.not.btdb_put(rtdb,rtdb_unita,mt_dbl,9,unita)) 2576 > call errquit('control_set_unita:writing unita',0,RTDB_ERR) 2577 end if 2578 value=btdb_parallel(.true.) 2579 return 2580 end 2581 2582 2583* *********************************** 2584* * * 2585* * control_unita_frozen * 2586* * * 2587* *********************************** 2588 real*8 function control_unita_frozen(i,j) 2589 implicit none 2590 integer i,j 2591 2592#include "control.fh" 2593 2594 control_unita_frozen = unita_frozen(i,j) 2595 return 2596 end 2597 2598* *********************************** 2599* * * 2600* * control_frozen * 2601* * * 2602* *********************************** 2603 logical function control_frozen() 2604 implicit none 2605 2606#include "control.fh" 2607 2608 control_frozen = frozen 2609 return 2610 end 2611 2612 2613* *********************************** 2614* * * 2615* * control_boundry * 2616* * * 2617* *********************************** 2618 character*12 function control_boundry() 2619 implicit none 2620 2621#include "control.fh" 2622 2623 control_boundry = boundry 2624 return 2625 end 2626 2627 2628c* *********************************** 2629c* * * 2630c* * control_pspnames * 2631c* * * 2632c* *********************************** 2633c character*20 function control_pspnames(i) 2634c implicit none 2635c integer i 2636c 2637c#include "control.fh" 2638c 2639c control_pspnames = pspnames(i) 2640c return 2641c end 2642c 2643c 2644c* *********************************** 2645c* * * 2646c* * control_pspstressnames * 2647c* * * 2648c* *********************************** 2649c character*20 function control_pspstressnames(i) 2650c implicit none 2651c integer i 2652c 2653c integer ind 2654c character*20 pspname 2655c character*20 control_pspnames 2656c external control_pspnames 2657c 2658c pspname = control_pspnames(i) 2659c ind = index(pspname,' ') -1 2660c pspname = pspname(1:ind)//'2' 2661c 2662c control_pspstressnames = pspname 2663c return 2664c end 2665 2666* *********************************** 2667* * * 2668* * control_npsp * 2669* * * 2670* *********************************** 2671 integer function control_npsp() 2672 implicit none 2673 2674#include "control.fh" 2675 2676 control_npsp = npsp 2677 return 2678 end 2679 2680 2681 2682* *********************************** 2683* * * 2684* * control_ecut * 2685* * * 2686* *********************************** 2687 real*8 function control_ecut() 2688 implicit none 2689 2690#include "control.fh" 2691 control_ecut = ecut 2692 return 2693 end 2694 2695 2696 2697* *********************************** 2698* * * 2699* * control_wcut * 2700* * * 2701* *********************************** 2702 real*8 function control_wcut() 2703 implicit none 2704 2705#include "control.fh" 2706 2707 control_wcut = wcut 2708 return 2709 end 2710 2711 2712* *********************************** 2713* * * 2714* * control_rcut * 2715* * * 2716* *********************************** 2717 real*8 function control_rcut() 2718 implicit none 2719 2720#include "control.fh" 2721 2722 control_rcut = rcut 2723 return 2724 end 2725 2726* *********************************** 2727* * * 2728* * control_ncut * 2729* * * 2730* *********************************** 2731 integer function control_ncut() 2732 implicit none 2733 2734#include "control.fh" 2735 2736 control_ncut = ncut 2737 return 2738 end 2739 2740 2741 2742 2743* *********************************** 2744* * * 2745* * control_output_psi * 2746* * * 2747* *********************************** 2748 character*50 function control_output_psi() 2749 implicit none 2750 2751#include "control.fh" 2752 2753 control_output_psi = output_wavefunction_filename 2754 return 2755 end 2756 2757* *********************************** 2758* * * 2759* * control_output_epsi * 2760* * * 2761* *********************************** 2762 character*50 function control_output_epsi() 2763 implicit none 2764 2765#include "control.fh" 2766 2767 control_output_epsi = output_ewavefunction_filename 2768 return 2769 end 2770 2771 2772* *********************************** 2773* * * 2774* * control_input_psi * 2775* * * 2776* *********************************** 2777 character*50 function control_input_psi() 2778 implicit none 2779 2780#include "control.fh" 2781 2782 control_input_psi = input_wavefunction_filename 2783 return 2784 end 2785 2786 2787* *********************************** 2788* * * 2789* * control_input_epsi * 2790* * * 2791* *********************************** 2792 character*50 function control_input_epsi() 2793 implicit none 2794 2795#include "control.fh" 2796 2797 control_input_epsi = input_ewavefunction_filename 2798 return 2799 end 2800 2801 2802* *********************************** 2803* * * 2804* * control_output_v_psi * 2805* * * 2806* *********************************** 2807 character*50 function control_output_v_psi() 2808 implicit none 2809 2810#include "control.fh" 2811 2812 control_output_v_psi = output_v_wavefunction_filename 2813 return 2814 end 2815 2816 2817* *********************************** 2818* * * 2819* * control_input_v_psi * 2820* * * 2821* *********************************** 2822 character*50 function control_input_v_psi() 2823 implicit none 2824 2825#include "control.fh" 2826 2827 control_input_v_psi = input_v_wavefunction_filename 2828 return 2829 end 2830 2831 2832* *********************************** 2833* * * 2834* * control_use_fractional_rho * 2835* * * 2836* *********************************** 2837 logical function control_use_fractional_rho() 2838 implicit none 2839 2840#include "bafdecls.fh" 2841#include "btdb.fh" 2842#include "errquit.fh" 2843#include "control.fh" 2844 2845* **** control_rtdb common block **** 2846 integer rtdb 2847 common / control_rtdb1 / rtdb 2848 2849 logical use_rho 2850 2851 if (.not.btdb_get(rtdb,"nwpw:use_fractional_rho", 2852 > mt_log,1,use_rho)) then 2853 use_rho = .true. 2854 end if 2855 control_use_fractional_rho = use_rho 2856 return 2857 end 2858 2859 2860* *********************************** 2861* * * 2862* * control_input_rho * 2863* * * 2864* *********************************** 2865 2866 character*50 function control_input_rho() 2867 implicit none 2868 2869#include "bafdecls.fh" 2870#include "btdb.fh" 2871#include "errquit.fh" 2872#include "control.fh" 2873 2874* **** control_rtdb common block **** 2875 integer rtdb 2876 common / control_rtdb1 / rtdb 2877 2878 character*50 filename 2879 2880 if (.not.btdb_cget(rtdb,'nwpw:input rho',1,filename)) then 2881 call util_file_prefix('rho',filename) 2882 end if 2883 control_input_rho = filename 2884 return 2885 end 2886 2887* *********************************** 2888* * * 2889* * control_output_rho * 2890* * * 2891* *********************************** 2892 2893 character*50 function control_output_rho() 2894 implicit none 2895 2896#include "bafdecls.fh" 2897#include "btdb.fh" 2898#include "errquit.fh" 2899#include "control.fh" 2900 2901* **** control_rtdb common block **** 2902 integer rtdb 2903 common / control_rtdb1 / rtdb 2904 2905 character*50 filename 2906 2907 if (.not.btdb_cget(rtdb,'nwpw:output rho',1,filename)) then 2908 call util_file_prefix('rho',filename) 2909 end if 2910 control_output_rho = filename 2911 return 2912 end 2913 2914 2915 2916 2917* *********************************** 2918* * * 2919* * control_xyz * 2920* * * 2921* *********************************** 2922 character*50 function control_xyz() 2923 implicit none 2924 2925 2926#include "control.fh" 2927 2928 control_xyz = xyz_filename 2929 return 2930 end 2931 2932* *********************************** 2933* * * 2934* * control_cell_name * 2935* * * 2936* *********************************** 2937 character*50 function control_cell_name() 2938 implicit none 2939 2940 2941#include "control.fh" 2942 2943 control_cell_name = cell_name 2944 return 2945 end 2946 2947 2948 2949* *********************************** 2950* * * 2951* * control_gga * 2952* * * 2953* *********************************** 2954 integer function control_gga() 2955 implicit none 2956 2957#include "control.fh" 2958 2959 control_gga = gga 2960 return 2961 end 2962 2963 2964* *********************************** 2965* * * 2966* * control_is_grimme2 * 2967* * * 2968* *********************************** 2969 logical function control_is_grimme2() 2970 implicit none 2971 2972#include "control.fh" 2973 2974 control_is_grimme2 = is_grimme2 2975 return 2976 end 2977 2978 2979* *********************************** 2980* * * 2981* * control_has_disp * 2982* * * 2983* *********************************** 2984 logical function control_has_disp() 2985 implicit none 2986 2987#include "control.fh" 2988 2989 control_has_disp = has_disp 2990 return 2991 end 2992 2993 2994* *********************************** 2995* * * 2996* * control_options_disp * 2997* * * 2998* *********************************** 2999 character*80 function control_options_disp() 3000 implicit none 3001 3002#include "control.fh" 3003 3004 control_options_disp = options_disp 3005 return 3006 end 3007 3008* *********************************** 3009* * * 3010* * control_has_vdw * 3011* * * 3012* *********************************** 3013 logical function control_has_vdw() 3014 implicit none 3015 3016#include "control.fh" 3017 3018 control_has_vdw = has_vdw 3019 return 3020 end 3021 3022* *********************************** 3023* * * 3024* * control_is_vdw2 * 3025* * * 3026* *********************************** 3027 logical function control_is_vdw2() 3028 implicit none 3029 3030#include "control.fh" 3031 3032 control_is_vdw2 = is_vdw2 3033 return 3034 end 3035 3036 3037* *********************************** 3038* * * 3039* * control_multiplicity * 3040* * * 3041* *********************************** 3042 integer function control_multiplicity() 3043 implicit none 3044 3045#include "control.fh" 3046 3047 control_multiplicity = multiplicity 3048 return 3049 end 3050 3051* *********************************** 3052* * * 3053* * control_multiplicity_set * 3054* * * 3055* *********************************** 3056 subroutine control_multiplicity_set(new_multiplicity) 3057 implicit none 3058 integer new_multiplicity 3059 3060#include "control.fh" 3061 3062 multiplicity = new_multiplicity 3063 return 3064 end 3065 3066* ***************************************** 3067* * * 3068* * control_check_charge_multiplicity * 3069* * * 3070* ***************************************** 3071 logical function control_check_charge_multiplicity() 3072 implicit none 3073 3074#include "control.fh" 3075 3076* *** local variables *** 3077 logical check 3078 real*8 icharge,tcharge,t 3079 integer mult,x,x_wf,nextra_orbs 3080 integer ispin_wf,ne_wf(2),x_f 3081 3082* ***** local functions **** 3083 logical psi_filefind 3084 external psi_filefind 3085 real*8 control_TotalCharge 3086 external control_TotalCharge 3087 real*8 ion_TotalCharge_qm 3088 external ion_TotalCharge_qm 3089 !integer control_frac_occ_extra_orbitals 3090 !external control_frac_occ_extra_orbitals 3091 integer control_fractional_orbitals 3092 external control_fractional_orbitals 3093 3094 3095 3096* **** check wavefunction file **** 3097 if (psi_filefind()) then 3098 3099* **** get mult and e-charge from wavefunction file **** 3100 call psi_get_ne(ispin_wf,ne_wf) 3101 3102c write(*,*) "ne_wf = ",ne_wf 3103c x_wf0 = ne_wf(1)+ne_wf(2) 3104c mult0 = ne_wf(1)-ne_wf(2) + 1 3105 !nextra_orbs = control_frac_occ_extra_orbitals() 3106 if (ispin_wf.eq.1) then 3107 ne_wf(1) = ne_wf(1) - control_fractional_orbitals(1) 3108 else 3109 ne_wf(1) = ne_wf(1) - control_fractional_orbitals(1) 3110 ne_wf(2) = ne_wf(2) - control_fractional_orbitals(2) 3111 end if 3112c write(*,*) "after ne_wf = ",ne_wf 3113 3114 x_wf = ne_wf(1)+ne_wf(2) 3115 mult = ne_wf(1)-ne_wf(2) + 1 3116 if (ispin_wf.eq.1) then 3117c x_wf0 = 2*x_wf0 3118c mult0 = 1 3119 x_wf = 2*x_wf 3120 mult = 1 3121 end if 3122 3123* **** get mult and e-charge from control **** 3124 tcharge = control_TotalCharge() 3125 icharge = ion_TotalCharge_qm() 3126 t = icharge - tcharge !** total number of electrons ** 3127 x = int(NINT(t)) 3128 3129 3130* **** reassign spin to agree with total number of electrons **** 3131 if ((mod(x,2).ne.0).and.(ispin.eq.1)) then !** odd number of electrons ** 3132 ispin = 2 3133 end if 3134 3135* **** reassign multiplicity to agree with total number of electrons **** 3136* *** odd number of electrons and mult odd *** 3137 if ((mod(x,2).ne.0) .and.(mod(multiplicity,2).ne.0)) then 3138 multiplicity = multiplicity - 1 3139 do while (multiplicity.gt.(x+1)) 3140 multiplicity = multiplicity - 2 3141 end do 3142 if (multiplicity.lt.1) multiplicity = 2 3143 end if 3144* *** even number of electrons and mult even *** 3145 if ((mod(x,2).eq.0) .and.(mod(multiplicity,2).eq.0)) then 3146 multiplicity = multiplicity - 1 3147 do while (multiplicity.gt.(x+1)) 3148 multiplicity = multiplicity - 2 3149 end do 3150 if (multiplicity.lt.1) multiplicity = 1 3151 end if 3152 3153 3154* **** compare multiplicity, charge, and ispin **** 3155c check = ( (((mult.eq.multiplicity).and.(x_wf.eq.x)) .or. 3156c > ((mult0.eq.multiplicity).and.(x_wf0.eq.x))).and. 3157c > (ispin.eq.ispin_wf)) 3158 check = ( (mult.eq.multiplicity).and.(x_wf.eq.x).and. 3159 > (ispin.eq.ispin_wf)) 3160 3161 if (ispin_wf.eq.3) then 3162 check=((ne_wf(1).eq.ne_wf(2)).and. 3163 > (x.eq.ne_wf(1) ) ) 3164 end if 3165* **** no wavefunction file *** 3166 else 3167 check = .false. 3168 end if 3169 3170 control_check_charge_multiplicity = check 3171 return 3172 end 3173 3174 3175 3176* ***************************************** 3177* * * 3178* * control_check_number_virtuals * 3179* * * 3180* ***************************************** 3181 logical function control_check_number_virtuals() 3182 implicit none 3183 3184#include "control.fh" 3185 3186* *** local variables *** 3187 logical check 3188 integer ispin_wf,ne_wf(2),ne(2) 3189 3190* ***** local functions **** 3191 logical epsi_filefind 3192 external epsi_filefind 3193 integer control_excited_ne 3194 external control_excited_ne 3195 3196* **** check wavefunction file **** 3197 if (epsi_filefind()) then 3198 3199* **** get mult and e-charge from wavefunction file **** 3200 call psi_get_ne_excited(ispin_wf,ne_wf) 3201 ne(1) = 0 3202 ne(2) = 0 3203 ne(1) = control_excited_ne(1) 3204 if (ispin.eq.2) ne(2) = control_excited_ne(2) 3205 3206 check = ((ne(1).eq.ne_wf(1)).and. 3207 > (ne(2).eq.ne_wf(2)).and. 3208 > (ispin.eq.ispin_wf)) 3209 3210* **** no wavefunction file *** 3211 else 3212 check = .false. 3213 end if 3214 3215 control_check_number_virtuals = check 3216 return 3217 end 3218 3219 3220 3221* *********************************** 3222* * * 3223* * control_ispin * 3224* * * 3225* *********************************** 3226 integer function control_ispin() 3227 implicit none 3228 3229#include "control.fh" 3230 3231 control_ispin = ispin 3232 return 3233 end 3234 3235* *********************************** 3236* * * 3237* * control_ispin_set * 3238* * * 3239* *********************************** 3240 subroutine control_ispin_set(new_ispin) 3241 implicit none 3242 integer new_ispin 3243 3244#include "control.fh" 3245 3246 ispin = new_ispin 3247 return 3248 end 3249 3250 3251* ******************************************* 3252* * * 3253* * control_gradient_iterations * 3254* * * 3255* ******************************************* 3256 subroutine control_gradient_iterations() 3257 implicit none 3258 3259#include "control.fh" 3260 3261 loop(1) = 1 3262 loop(2) = 1 3263 3264 return 3265 end 3266 3267* *********************************** 3268* * * 3269* * control_version * 3270* * * 3271* *********************************** 3272 integer function control_version() 3273 implicit none 3274 3275#include "inp.fh" 3276#include "control.fh" 3277 3278* **** local variables **** 3279 integer l,version 3280 3281 l =index(boundry,' ') - 1 3282 3283 version = 3 3284 if (inp_compare(.false.,boundry(1:l),'periodic')) version=3 3285 if (inp_compare(.false.,boundry(1:l),'aperiodic')) version=4 3286 3287 control_version = version 3288 return 3289 end 3290 3291 3292* ************************ 3293* * * 3294* * control_Nose * 3295* * * 3296* ************************ 3297 logical function control_Nose() 3298 implicit none 3299 3300* **** control_nose common block **** 3301 logical nose 3302 real*8 Pe,Te,Pr,Tr 3303 common / control_nblock / Pe,Te,Pr,Tr,nose 3304 3305 3306 control_Nose = nose 3307 return 3308 end 3309 3310 3311* **************************** 3312* * * 3313* * control_Nose_Pe * 3314* * * 3315* **************************** 3316 real*8 function control_Nose_Pe() 3317 implicit none 3318 3319* **** control_nose common block **** 3320 logical nose 3321 real*8 Pe,Te,Pr,Tr 3322 common / control_nblock / Pe,Te,Pr,Tr,nose 3323 3324 3325 control_Nose_Pe = Pe 3326 return 3327 end 3328 3329 3330* **************************** 3331* * * 3332* * control_Nose_Te * 3333* * * 3334* **************************** 3335 real*8 function control_Nose_Te() 3336 implicit none 3337 3338* **** control_nose common block **** 3339 logical nose 3340 real*8 Pe,Te,Pr,Tr 3341 common / control_nblock / Pe,Te,Pr,Tr,nose 3342 3343 3344 control_Nose_Te = Te 3345 return 3346 end 3347 3348 3349* **************************** 3350* * * 3351* * control_Nose_Pr * 3352* * * 3353* **************************** 3354 real*8 function control_Nose_Pr() 3355 implicit none 3356 3357* **** control_nose common block **** 3358 logical nose 3359 real*8 Pe,Te,Pr,Tr 3360 common / control_nblock / Pe,Te,Pr,Tr,nose 3361 3362 3363 control_Nose_Pr = Pr 3364 return 3365 end 3366 3367 3368* **************************** 3369* * * 3370* * control_Nose_Tr * 3371* * * 3372* **************************** 3373 real*8 function control_Nose_Tr() 3374 implicit none 3375 3376* **** control_nose common block **** 3377 logical nose 3378 real*8 Pe,Te,Pr,Tr 3379 common / control_nblock / Pe,Te,Pr,Tr,nose 3380 3381 3382 control_Nose_Tr = Tr 3383 return 3384 end 3385 3386 3387* **************************** 3388* * * 3389* * control_Mulliken * 3390* * * 3391* **************************** 3392 logical function control_Mulliken() 3393 implicit none 3394 3395#include "bafdecls.fh" 3396#include "btdb.fh" 3397#include "errquit.fh" 3398#include "control.fh" 3399 3400* **** control_rtdb common block **** 3401 integer rtdb 3402 common / control_rtdb1 / rtdb 3403 3404* ***** local variables **** 3405 logical value 3406 3407 value = .false. 3408 if (code.eq.1) then 3409 if (.not.btdb_get(rtdb,'cpsd:mulliken',mt_log,1,value)) 3410 > value = .false. 3411 end if 3412 if ((code.eq.2).or.(code.eq.7)) then 3413 if (.not.btdb_get(rtdb,'cpmd:mulliken',mt_log,1,value)) 3414 > value = .false. 3415 end if 3416 3417 3418 if ((code.eq.3).or.(code.eq.11)) then 3419 if (.not.btdb_get(rtdb,'cgsd:mulliken',mt_log,1,value)) 3420 > value = .false. 3421 end if 3422 3423 if ((code.eq.5).or.(code.eq.11)) then 3424 if (.not.btdb_get(rtdb,'band:mulliken',mt_log,1,value)) 3425 > value = .false. 3426 end if 3427 3428 3429 control_Mulliken = value 3430 return 3431 end 3432 3433 3434* ***************************** 3435* * * 3436* * control_allow_translation * 3437* * * 3438* ***************************** 3439 logical function control_allow_translation() 3440 implicit none 3441 3442#include "bafdecls.fh" 3443#include "btdb.fh" 3444#include "inp.fh" 3445#include "util.fh" 3446 3447 logical value 3448 character*50 operation 3449 3450 logical control_qmmm 3451 external control_qmmm 3452 3453* **** control_rtdb common block **** 3454 integer rtdb 3455 common / control_rtdb1 / rtdb 3456 3457 if (.not.btdb_get(rtdb,'cgsd:allow_translation', 3458 > mt_log,1,value)) 3459 > value = .true. 3460 3461* *** read the current operation **** 3462 if (.not. btdb_cget(rtdb, 'task:operation', 1, operation)) 3463 $ operation = ' ' 3464 3465* *** allow translation of operation == freq||hessian *** 3466 if (inp_compare(.false.,'freq',operation)) value = .true. 3467 if (inp_compare(.false.,'hessian',operation)) value = .true. 3468 3469* *** allow translation for QMMM *** 3470 if (control_qmmm()) value = .true. 3471 3472 3473 control_allow_translation = value 3474 return 3475 end 3476 3477 3478 3479 3480* ***************************** 3481* * * 3482* * control_qmmm * 3483* * * 3484* ***************************** 3485 logical function control_qmmm() 3486 implicit none 3487 3488#include "bafdecls.fh" 3489#include "btdb.fh" 3490#include "util.fh" 3491 3492 logical task_qmmm 3493 3494* **** control_rtdb common block **** 3495 integer rtdb 3496 common / control_rtdb1 / rtdb 3497 3498 if( .not. btdb_get(rtdb,'task:QMMM',mt_log,1,task_qmmm)) 3499 > task_qmmm = .false. 3500 3501 control_qmmm = task_qmmm 3502 return 3503 end 3504 3505 3506* ***************************** 3507* * * 3508* * control_makehmass2 * 3509* * * 3510* ***************************** 3511 logical function control_makehmass2() 3512 implicit none 3513 3514#include "bafdecls.fh" 3515#include "btdb.fh" 3516#include "util.fh" 3517 3518 logical makehmass2 3519 3520* **** control_rtdb common block **** 3521 integer rtdb 3522 common / control_rtdb1 / rtdb 3523 3524 if( .not. btdb_get(rtdb,'nwpw:makehmass2',mt_log,1,makehmass2)) 3525 > makehmass2 = .true. 3526 3527 control_makehmass2 = makehmass2 3528 return 3529 end 3530 3531 3532 3533* ***************************** 3534* * * 3535* * control_MP2 * 3536* * * 3537* ***************************** 3538 logical function control_MP2() 3539 implicit none 3540 3541#include "bafdecls.fh" 3542#include "btdb.fh" 3543#include "util.fh" 3544 3545 logical mp2 3546 3547* **** control_rtdb common block **** 3548 integer rtdb 3549 common / control_rtdb1 / rtdb 3550 3551 if( .not. btdb_get(rtdb,'nwpw:MP2',mt_log,1,mp2)) 3552 > mp2 = .false. 3553 3554 control_MP2 = mp2 3555 return 3556 end 3557 3558 3559* ***************************** 3560* * * 3561* * control_2qintegrals * 3562* * * 3563* ***************************** 3564 logical function control_2qintegrals() 3565 implicit none 3566 3567#include "bafdecls.fh" 3568#include "btdb.fh" 3569#include "util.fh" 3570 3571 logical q2 3572 3573* **** control_rtdb common block **** 3574 integer rtdb 3575 common / control_rtdb1 / rtdb 3576 3577 if( .not. btdb_get(rtdb,'nwpw:2qintegrals',mt_log,1,q2)) 3578 > q2 = .false. 3579 3580 control_2qintegrals = q2 3581 return 3582 end 3583 3584 3585* **************************** 3586* * * 3587* * control_num_kvectors * 3588* * * 3589* **************************** 3590 integer function control_num_kvectors() 3591 implicit none 3592 3593#include "bafdecls.fh" 3594#include "btdb.fh" 3595#include "errquit.fh" 3596 3597* **** control_rtdb common block **** 3598 integer rtdb 3599 common / control_rtdb1 / rtdb 3600 3601* **** local variables **** 3602 logical value 3603 character*50 zone_name 3604 character*50 rtdb_name 3605 integer num_kvectors,l 3606 3607 if (.not.btdb_cget(rtdb,'band:zone_name',1,zone_name)) 3608 > zone_name = 'zone_default' 3609 3610 l = index(zone_name,' ') -1 3611 rtdb_name = zone_name(1:l)//':number_kvectors' 3612 if (.not.btdb_get(rtdb,rtdb_name,mt_int,1,num_kvectors)) 3613 > num_kvectors = 1 3614 3615 control_num_kvectors = num_kvectors 3616 return 3617 end 3618 3619* ********************************** 3620* * * 3621* * control_num_kvectors_print * 3622* * * 3623* ********************************** 3624 integer function control_num_kvectors_print() 3625 implicit none 3626 3627#include "bafdecls.fh" 3628#include "btdb.fh" 3629#include "errquit.fh" 3630 3631* **** control_rtdb common block **** 3632 integer rtdb 3633 common / control_rtdb1 / rtdb 3634 3635* **** local variables **** 3636 logical value 3637 character*50 zone_name 3638 character*50 rtdb_name 3639 integer num_kvectors_print,l 3640 3641 if (.not.btdb_cget(rtdb,'band:zone_name',1,zone_name)) 3642 > zone_name = 'zone_default' 3643 3644 l = index(zone_name,' ') -1 3645 rtdb_name = zone_name(1:l)//':number_kvectors_print' 3646 if (.not.btdb_get(rtdb,rtdb_name,mt_int,1,num_kvectors_print)) 3647 > num_kvectors_print = 0 3648 3649 control_num_kvectors_print = num_kvectors_print 3650 return 3651 end 3652 3653 3654* **************************** 3655* * * 3656* * control_ksvector * 3657* * * 3658* **************************** 3659 subroutine control_ksvector(i,ks) 3660 implicit none 3661 integer i 3662 real*8 ks(4) 3663 3664#include "bafdecls.fh" 3665#include "btdb.fh" 3666#include "errquit.fh" 3667 3668* **** control_rtdb common block **** 3669 integer rtdb 3670 common / control_rtdb1 / rtdb 3671 3672* **** local variables **** 3673 character*50 zone_name 3674 character*50 rtdb_name 3675 integer num_kvectors,l 3676 integer kvs(2) 3677 3678* **** external functions **** 3679 integer control_num_kvectors 3680 external control_num_kvectors 3681 3682 num_kvectors = control_num_kvectors() 3683 3684 if (.not.BA_push_get(mt_dbl,(4*num_kvectors),'kvs',kvs(2),kvs(1))) 3685 > call errquit('control_ksvector: out of stack', 0,MA_ERR) 3686 3687 if (.not.btdb_cget(rtdb,'band:zone_name',1,zone_name)) 3688 > zone_name = 'zone_default' 3689 3690 3691 l = index(zone_name,' ') -1 3692 rtdb_name = zone_name(1:l)//':kvectors' 3693 if (.not.btdb_get(rtdb,rtdb_name,mt_dbl, 3694 > (4*num_kvectors),dbl_mb(kvs(1)))) then 3695 call dcopy(4*num_kvectors,0.0d0,0,dbl_mb(kvs(1)),1) 3696 end if 3697 3698 ks(1) = dbl_mb(kvs(1)+4*(i-1)) 3699 ks(2) = dbl_mb(kvs(1)+4*(i-1)+1) 3700 ks(3) = dbl_mb(kvs(1)+4*(i-1)+2) 3701 ks(4) = dbl_mb(kvs(1)+4*(i-1)+3) 3702 3703 if (.not.BA_pop_stack(kvs(2))) 3704 > call errquit('control_ksvector: failed to free stack',0,MA_ERR) 3705 3706 return 3707 end 3708 3709* **************************** 3710* * * 3711* * control_kvector * 3712* * * 3713* **************************** 3714 subroutine control_kvector(i,kv) 3715 implicit none 3716 integer i 3717 real*8 kv(3) 3718 3719* **** local variables **** 3720 real*8 ks(4) 3721 3722* **** external functions **** 3723 real*8 lattice_unitg 3724 external lattice_unitg 3725 3726 call control_ksvector(i,ks) 3727 3728 kv(1) = ks(1)*lattice_unitg(1,1) 3729 > + ks(2)*lattice_unitg(1,2) 3730 > + ks(3)*lattice_unitg(1,3) 3731 kv(2) = ks(1)*lattice_unitg(2,1) 3732 > + ks(2)*lattice_unitg(2,2) 3733 > + ks(3)*lattice_unitg(2,3) 3734 kv(3) = ks(1)*lattice_unitg(3,1) 3735 > + ks(2)*lattice_unitg(3,2) 3736 > + ks(3)*lattice_unitg(3,3) 3737 3738 return 3739 end 3740 3741* ********************************** 3742* * * 3743* * control_monkhorst_pack_grid * 3744* * * 3745* ********************************** 3746 subroutine control_monkhorst_pack_grid(grid) 3747 implicit none 3748 integer grid(3) 3749 3750#include "bafdecls.fh" 3751#include "btdb.fh" 3752#include "errquit.fh" 3753 3754* **** control_rtdb common block **** 3755 integer rtdb 3756 common / control_rtdb1 / rtdb 3757 3758* **** local variables **** 3759 character*50 zone_name 3760 character*50 rtdb_name 3761 integer l 3762 3763 if (.not.btdb_cget(rtdb,'band:zone_name',1,zone_name)) 3764 > zone_name = 'zone_default' 3765 3766 l = index(zone_name,' ') -1 3767 rtdb_name = zone_name(1:l)//':monkhorst-pack' 3768 if (.not.btdb_get(rtdb,rtdb_name,mt_int,3,grid)) then 3769 grid(1) = 0 3770 grid(2) = 0 3771 grid(3) = 0 3772 end if 3773 3774 return 3775 end 3776 3777 3778* ********************************** 3779* * * 3780* * control_ksvector_index * 3781* * * 3782* ********************************** 3783 integer function control_ksvector_index(ks) 3784 implicit none 3785 real*8 ks(*) 3786 3787#include "bafdecls.fh" 3788#include "btdb.fh" 3789#include "errquit.fh" 3790 3791* **** control_rtdb common block **** 3792 integer rtdb 3793 common / control_rtdb1 / rtdb 3794 3795* **** local variables **** 3796 character*50 zone_name 3797 character*50 rtdb_name 3798 integer num_kvectors,i,i1,i2 3799 integer kvs(2) 3800 real*8 d1,d2 3801 3802* **** external functions **** 3803 integer control_num_kvectors 3804 external control_num_kvectors 3805 3806 num_kvectors = control_num_kvectors() 3807 3808 if (.not.BA_push_get(mt_dbl,(4*num_kvectors),'kvs',kvs(2),kvs(1))) 3809 > call errquit('control_ksvector: out of stack', 0,MA_ERR) 3810 3811 if (.not.btdb_cget(rtdb,'band:zone_name',1,zone_name)) 3812 > zone_name = 'zone_default' 3813 3814 i = index(zone_name,' ') -1 3815 rtdb_name = zone_name(1:i)//':kvectors' 3816 if (.not.btdb_get(rtdb,rtdb_name,mt_dbl, 3817 > (4*num_kvectors),dbl_mb(kvs(1)))) then 3818 call dcopy(4*num_kvectors,0.0d0,0,dbl_mb(kvs(1)),1) 3819 end if 3820 3821 i1 = -1 3822 i2 = -1 3823 do i=1,num_kvectors 3824 d1 = dsqrt((dbl_mb(kvs(1)+4*(i-1)) - ks(1))**2 3825 > + (dbl_mb(kvs(1)+4*(i-1)+1) - ks(2))**2 3826 > + (dbl_mb(kvs(1)+4*(i-1)+2) - ks(3))**2) 3827 if (d1.lt.1.0e-6) i1 = i 3828 d2 = dsqrt((dbl_mb(kvs(1)+4*(i-1)) + ks(1))**2 3829 > + (dbl_mb(kvs(1)+4*(i-1)+1) + ks(2))**2 3830 > + (dbl_mb(kvs(1)+4*(i-1)+2) + ks(3))**2) 3831 if (d2.lt.1.0e-6) i2 = i 3832 end do 3833 3834 if (.not.BA_pop_stack(kvs(2))) 3835 > call errquit('control_ksvector: failed to free stack',0,MA_ERR) 3836 3837 if (i1.eq.-1) i1 = i2 3838 control_ksvector_index = i1 3839 return 3840 end 3841 3842 3843 3844 3845* ***************************** 3846* * * 3847* * control_TotalCharge * 3848* * * 3849* ***************************** 3850 real*8 function control_TotalCharge() 3851 implicit none 3852 3853#include "bafdecls.fh" 3854#include "btdb.fh" 3855 3856* **** control_rtdb common block **** 3857 integer rtdb 3858 common / control_rtdb1 / rtdb 3859 3860 double precision charge 3861 3862 charge = 0.0d0 3863 if (.not.btdb_get(rtdb,'charge',mt_dbl,1,charge)) then 3864 charge = 0.0d0 3865 end if 3866 3867 control_TotalCharge = charge 3868 return 3869 end 3870 3871 3872 3873* ************************************** 3874* * * 3875* * control_frac_occ_extra_orbitals * 3876* * * 3877* ************************************** 3878 integer function control_frac_occ_extra_orbitals() 3879 implicit none 3880 3881#include "bafdecls.fh" 3882#include "btdb.fh" 3883 3884* **** control_rtdb common block **** 3885 integer rtdb 3886 common / control_rtdb1 / rtdb 3887 3888 integer norbs 3889 3890 norbs = 0 3891 if (.not. 3892 > btdb_get(rtdb,'nwpw:frac_occ:extra_orbitals', 3893 > mt_int,1,norbs)) then 3894 norbs = 0 3895 end if 3896 3897 control_frac_occ_extra_orbitals = norbs 3898 return 3899 end 3900 3901 3902 3903 3904* ***************************** 3905* * * 3906* * control_rtdb * 3907* * * 3908* ***************************** 3909 integer function control_rtdb() 3910 implicit none 3911 3912* **** control_rtdb common block **** 3913 integer trtdb 3914 common / control_rtdb1 / trtdb 3915 3916 control_rtdb = trtdb 3917 return 3918 end 3919 3920* ***************************** 3921* * * 3922* * control_minimizer * 3923* * * 3924* ***************************** 3925 integer function control_minimizer() 3926 implicit none 3927 3928#include "control.fh" 3929 3930 control_minimizer = minimizer 3931 return 3932 end 3933 3934 3935 3936* ***************************** 3937* * * 3938* * control_lmbfgs_size * 3939* * * 3940* ***************************** 3941 integer function control_lmbfgs_size() 3942 implicit none 3943 3944#include "bafdecls.fh" 3945#include "btdb.fh" 3946 3947* **** control_rtdb common block **** 3948 integer rtdb 3949 common / control_rtdb1 / rtdb 3950 3951 integer lmbfgs_size 3952 3953 if (.not.btdb_get(rtdb,'nwpw:lmbfgs_size',mt_int,1,lmbfgs_size)) 3954 > then 3955 lmbfgs_size = 1 3956 end if 3957 3958 control_lmbfgs_size = lmbfgs_size 3959 return 3960 end 3961 3962* ***************************** 3963* * * 3964* * control_precondition * 3965* * * 3966* ***************************** 3967 logical function control_precondition() 3968 implicit none 3969 3970#include "bafdecls.fh" 3971#include "btdb.fh" 3972 3973* **** control_rtdb common block **** 3974 integer rtdb 3975 common / control_rtdb1 / rtdb 3976 3977 logical precondition 3978 3979 if (.not.btdb_get(rtdb,'nwpw:precondition', 3980 > mt_log,1,precondition)) 3981 > then 3982 precondition = .false. 3983 end if 3984 3985 control_precondition = precondition 3986 return 3987 end 3988 3989* ***************************** 3990* * * 3991* * control_lmbfgs_ondisk * 3992* * * 3993* ***************************** 3994 logical function control_lmbfgs_ondisk() 3995 implicit none 3996 3997#include "bafdecls.fh" 3998#include "btdb.fh" 3999 4000* **** control_rtdb common block **** 4001 integer rtdb 4002 common / control_rtdb1 / rtdb 4003 4004 logical ondisk 4005 4006 if (.not.btdb_get(rtdb,'nwpw:lmbfgs_ondisk', 4007 > mt_log,1,ondisk)) 4008 > then 4009 ondisk = .false. 4010 end if 4011 4012 control_lmbfgs_ondisk = ondisk 4013 return 4014 end 4015 4016* ***************************** 4017* * * 4018* * control_pspparameters * 4019* * * 4020* ***************************** 4021 subroutine control_pspparameters(atom,lmax,locp,rlocal) 4022 implicit none 4023 character*4 atom 4024 integer lmax,locp 4025 real*8 rlocal 4026 4027#include "bafdecls.fh" 4028#include "btdb.fh" 4029 4030* **** control_rtdb common block **** 4031 integer rtdb 4032 common / control_rtdb1 / rtdb 4033 4034 integer l,k 4035 character*5 element 4036 character*50 rtdb_name 4037 4038 element = ' ' 4039 element = atom 4040 l = index(element,' ') - 1 4041 rtdb_name = ' ' 4042 rtdb_name = element(1:l)//':lmax' 4043 k = index(rtdb_name,' ') - 1 4044 if (.not.btdb_get(rtdb,rtdb_name(1:k),mt_int,1,lmax)) 4045 > lmax = -1 4046 4047 rtdb_name = ' ' 4048 rtdb_name = element(1:l)//':locp' 4049 k = index(rtdb_name,' ') - 1 4050 if (.not.btdb_get(rtdb,rtdb_name(1:k),mt_int,1,locp)) 4051 > locp = -1 4052 4053 rtdb_name = ' ' 4054 rtdb_name = element(1:l)//':rlocal' 4055 k = index(rtdb_name,' ') - 1 4056 if (.not.btdb_get(rtdb,rtdb_name(1:k),mt_dbl,1,rlocal)) 4057 > rlocal = 1.0d0 4058 4059 return 4060 end 4061 4062* ********************************** 4063* * * 4064* * control_mullikenparameters * 4065* * * 4066* ********************************** 4067 subroutine control_mullikenparameters(atom,rcut,lmbda) 4068 implicit none 4069 character*4 atom 4070 real*8 rcut,lmbda 4071 4072#include "bafdecls.fh" 4073#include "btdb.fh" 4074 4075* **** control_rtdb common block **** 4076 integer rtdb 4077 common / control_rtdb1 / rtdb 4078 4079 logical value 4080 integer l,k 4081 character*5 element 4082 character*50 rtdb_name 4083 4084 element = ' ' 4085 element = atom 4086 l = index(element,' ') - 1 4087 rtdb_name = ' ' 4088 rtdb_name = element(1:l)//':mulliken:rcut' 4089 k = index(rtdb_name,' ') - 1 4090 if (.not.btdb_get(rtdb,rtdb_name(1:k),mt_dbl,1,rcut)) 4091 > rcut = 1.0d0 4092 4093 rtdb_name = ' ' 4094 rtdb_name = element(1:l)//':mulliken:lmbda' 4095 k = index(rtdb_name,' ') - 1 4096 if (.not.btdb_get(rtdb,rtdb_name(1:k),mt_dbl,1,lmbda)) 4097 > lmbda = 0.0d0 4098 4099 return 4100 end 4101 4102 4103 4104 4105* *************************** 4106* * * 4107* * control_Ep * 4108* * * 4109* *************************** 4110 real*8 function control_Ep() 4111 implicit none 4112 4113#include "bafdecls.fh" 4114#include "btdb.fh" 4115 4116* **** control_rtdb common block **** 4117 integer rtdb 4118 common / control_rtdb1 / rtdb 4119 4120 real*8 Ep 4121 4122 if (.not.btdb_get(rtdb,'nwpw:Eprecondition', 4123 > mt_dbl,1,Ep)) 4124 > then 4125 Ep = 20.0d0 4126 end if 4127 4128 control_Ep = Ep 4129 return 4130 end 4131 4132 4133 4134* *************************** 4135* * * 4136* * control_Sp * 4137* * * 4138* *************************** 4139 real*8 function control_Sp() 4140 implicit none 4141 4142#include "bafdecls.fh" 4143#include "btdb.fh" 4144 4145* **** control_rtdb common block **** 4146 integer rtdb 4147 common / control_rtdb1 / rtdb 4148 4149 real*8 Sp 4150 4151 if (.not.btdb_get(rtdb,'nwpw:Sprecondition', 4152 > mt_dbl,1,Sp)) 4153 > then 4154 Sp = 200.0d0 4155 end if 4156 4157 control_Sp = Sp 4158 return 4159 end 4160 4161 4162* *************************** 4163* * * 4164* * control_SA * 4165* * * 4166* *************************** 4167 logical function control_SA() 4168 implicit none 4169 4170#include "control.fh" 4171 4172 control_SA=SA 4173 return 4174 end 4175 4176* *************************** 4177* * * 4178* * control_SA_decay * 4179* * * 4180* *************************** 4181 real*8 function control_SA_decay(choice) 4182 implicit none 4183 integer choice 4184 4185#include "control.fh" 4186 4187 control_SA_decay = sa_decay(choice) 4188 return 4189 end 4190 4191* ********************************** 4192* * * 4193* * control_dipole_motion * 4194* * * 4195* ********************************** 4196 logical function control_dipole_motion() 4197 implicit none 4198 4199#include "control.fh" 4200 4201 control_dipole_motion=dipole_motion 4202 return 4203 end 4204 4205* *************************** 4206* * * 4207* * control_Fei * 4208* * * 4209* *************************** 4210 logical function control_Fei() 4211 implicit none 4212 4213#include "control.fh" 4214 4215 control_Fei=fei 4216 return 4217 end 4218 4219 4220* *************************** 4221* * * 4222* * control_Fei_quench * 4223* * * 4224* *************************** 4225 logical function control_Fei_quench() 4226 implicit none 4227 4228#include "control.fh" 4229 4230 control_Fei_quench=fei_quench 4231 return 4232 end 4233 4234 4235* *************************** 4236* * * 4237* * control_gram_schmidt * 4238* * * 4239* *************************** 4240 logical function control_gram_schmidt() 4241 implicit none 4242 4243#include "control.fh" 4244 4245 control_gram_schmidt=gram_schmidt 4246 return 4247 end 4248 4249 4250* *************************** 4251* * * 4252* * control_balance * 4253* * * 4254* *************************** 4255 logical function control_balance() 4256 implicit none 4257 4258#include "control.fh" 4259 4260 control_balance=balance 4261 return 4262 end 4263 4264 4265 4266* ***************************** 4267* * * 4268* * control_excited_ne * 4269* * * 4270* ***************************** 4271 integer function control_excited_ne(ii) 4272 implicit none 4273 integer ii 4274 4275#include "bafdecls.fh" 4276#include "btdb.fh" 4277 4278* **** control_rtdb common block **** 4279 integer rtdb 4280 common / control_rtdb1 / rtdb 4281 4282 integer ne(2) 4283 4284 if (.not.btdb_get(rtdb,'nwpw:excited_ne',mt_int,2,ne)) then 4285 ne(1) = 0 4286 ne(2) = 0 4287 end if 4288 4289 control_excited_ne = ne(ii) 4290 return 4291 end 4292 4293 4294 4295* ************************************* 4296* * * 4297* * control_fractional_orbitals * 4298* * * 4299* ************************************* 4300 integer function control_fractional_orbitals(ii) 4301 implicit none 4302 integer ii 4303#include "control.fh" 4304 control_fractional_orbitals = frac_ne(ii) 4305 return 4306 end 4307* ************************************* 4308* * * 4309* * control_fractional_kT * 4310* * * 4311* ************************************* 4312 real*8 function control_fractional_kT() 4313 implicit none 4314#include "control.fh" 4315* *** local variables and parameters **** 4316 double precision kb 4317 parameter (kb=3.16679d-6) 4318 control_fractional_kT = kb*frac_temperature 4319 return 4320 end 4321* ************************************* 4322* * * 4323* * control_fractional_temperature * 4324* * * 4325* ************************************* 4326 real*8 function control_fractional_temperature() 4327 implicit none 4328#include "control.fh" 4329 control_fractional_temperature = frac_temperature 4330 return 4331 end 4332* ************************************* 4333* * * 4334* * control_fractional_smeartype * 4335* * * 4336* ************************************* 4337 integer function control_fractional_smeartype() 4338 implicit none 4339#include "control.fh" 4340 control_fractional_smeartype = frac_smeartype 4341 return 4342 end 4343* ************************************* 4344* * * 4345* * control_fractional * 4346* * * 4347* ************************************* 4348 logical function control_fractional() 4349 implicit none 4350#include "control.fh" 4351 control_fractional = fractional 4352 return 4353 end 4354 4355* ************************************* 4356* * * 4357* * control_ortho * 4358* * * 4359* ************************************* 4360 logical function control_ortho() 4361 implicit none 4362 4363#include "bafdecls.fh" 4364#include "btdb.fh" 4365 4366* **** control_rtdb common block **** 4367 integer rtdb 4368 common / control_rtdb1 / rtdb 4369 4370 logical ortho 4371 4372 if (.not. 4373 > btdb_get(rtdb,'nwpw:ortho_initialize',mt_log,1,ortho)) then 4374 ortho = .true. 4375 end if 4376 4377 control_ortho = ortho 4378 return 4379 end 4380 4381 4382 4383 4384* ***************************** 4385* * * 4386* * control_mapping * 4387* * * 4388* ***************************** 4389 integer function control_mapping() 4390 implicit none 4391 4392#include "control.fh" 4393 4394 control_mapping = mapping 4395 return 4396 end 4397 4398 4399* ***************************** 4400* * * 4401* * control_np_dimensions * 4402* * * 4403* ***************************** 4404 integer function control_np_dimensions(i) 4405 implicit none 4406 integer i 4407 4408#include "control.fh" 4409 4410 control_np_dimensions = np_dimensions(i) 4411 return 4412 end 4413 4414* ***************************** 4415* * * 4416* * control_np_orbital * 4417* * * 4418* ***************************** 4419 integer function control_np_orbital() 4420 implicit none 4421 4422#include "control.fh" 4423 4424 control_np_orbital = np_dimensions(2) 4425 return 4426 end 4427 4428 4429 4430* ***************************** 4431* * * 4432* * control_mapping1d * 4433* * * 4434* ***************************** 4435 integer function control_mapping1d() 4436 implicit none 4437 4438#include "control.fh" 4439 4440 control_mapping1d = mapping1d 4441 return 4442 end 4443 4444* ***************************** 4445* * * 4446* * control_parallel_io * 4447* * * 4448* ***************************** 4449 logical function control_parallel_io() 4450 implicit none 4451 4452#include "control.fh" 4453 4454 control_parallel_io = pio 4455 return 4456 end 4457 4458 4459 4460* ***************************** 4461* * * 4462* * control_oep * 4463* * * 4464* ***************************** 4465 logical function control_oep() 4466 implicit none 4467 4468#include "bafdecls.fh" 4469#include "btdb.fh" 4470 4471* **** control_rtdb common block **** 4472 integer rtdb 4473 common / control_rtdb1 / rtdb 4474 4475 logical oep 4476 4477 if (.not.btdb_get(rtdb,'nwpw:oep',mt_log,1,oep)) then 4478 oep = .false. 4479 end if 4480 4481 control_oep = oep 4482 return 4483 end 4484 4485 4486 4487* ***************************** 4488* * * 4489* * control_new_vpsi * 4490* * * 4491* ***************************** 4492 logical function control_new_vpsi() 4493 implicit none 4494 4495#include "bafdecls.fh" 4496#include "btdb.fh" 4497 4498* **** control_rtdb common block **** 4499 integer rtdb 4500 common / control_rtdb1 / rtdb 4501 4502 logical new_vpsi 4503 4504 if (.not.btdb_get(rtdb,'nwpw:new_vmovecs',mt_log,1,new_vpsi)) then 4505 new_vpsi = .false. 4506 end if 4507 4508 control_new_vpsi = new_vpsi 4509 return 4510 end 4511 4512 4513* ***************************** 4514* * * 4515* * control_COM_shift * 4516* * * 4517* ***************************** 4518 logical function control_COM_shift() 4519 implicit none 4520 4521#include "bafdecls.fh" 4522#include "btdb.fh" 4523 4524* **** control_rtdb common block **** 4525 integer rtdb 4526 common / control_rtdb1 / rtdb 4527 4528 logical com_shift 4529 4530 if (.not.btdb_get(rtdb,'nwpw:com_shift',mt_log,1,com_shift)) then 4531 com_shift = .true. 4532 end if 4533 4534 control_COM_shift = com_shift 4535 return 4536 end 4537 4538 4539 4540* ***************************** 4541* * * 4542* * control_DOS * 4543* * * 4544* ***************************** 4545 logical function control_DOS() 4546 implicit none 4547 4548#include "bafdecls.fh" 4549#include "btdb.fh" 4550 4551* **** control_rtdb common block **** 4552 integer rtdb 4553 common / control_rtdb1 / rtdb 4554 4555 logical dos 4556 real*8 alpha 4557 4558 dos = .false. 4559 if (btdb_get(rtdb,'dos:alpha',mt_dbl,1,alpha)) dos = .true. 4560 4561 control_DOS = dos 4562 return 4563 end 4564 4565 4566* ***************************** 4567* * * 4568* * control_psi_tmp * 4569* * * 4570* ***************************** 4571 logical function control_psi_tmp() 4572 implicit none 4573 4574#include "bafdecls.fh" 4575#include "btdb.fh" 4576 4577* **** control_rtdb common block **** 4578 integer rtdb 4579 common / control_rtdb1 / rtdb 4580 4581 logical psitmp 4582 4583 psitmp = .false. 4584 if (.not.btdb_get(rtdb,'nwpw:psi_tmp',mt_log,1,psitmp)) 4585 > psitmp = .false. 4586 4587 control_psi_tmp = psitmp 4588 return 4589 end 4590 4591* ***************************** 4592* * * 4593* * control_mulliken_kawai * 4594* * * 4595* ***************************** 4596 logical function control_mulliken_kawai() 4597 implicit none 4598 4599#include "bafdecls.fh" 4600#include "btdb.fh" 4601 4602* **** control_rtdb common block **** 4603 integer rtdb 4604 common / control_rtdb1 / rtdb 4605 4606 logical value 4607 4608 value = .false. 4609 if (.not.btdb_get(rtdb,'nwpw:mulliken_kawai',mt_log,1,value)) 4610 > value = .false. 4611 4612 control_mulliken_kawai = value 4613 return 4614 end 4615 4616 4617 4618* ***************************** 4619* * * 4620* * control_zero_forces * 4621* * * 4622* ***************************** 4623 logical function control_zero_forces() 4624 implicit none 4625 4626#include "bafdecls.fh" 4627#include "btdb.fh" 4628 4629* **** control_rtdb common block **** 4630 integer rtdb 4631 common / control_rtdb1 / rtdb 4632 4633 logical value 4634 4635 value = .false. 4636 if (.not.btdb_get(rtdb,'nwpw:zero_forces',mt_log,1,value)) 4637 > value = .false. 4638 4639 control_zero_forces = value 4640 return 4641 end 4642 4643 4644 4645* ******************************** 4646* * * 4647* * control_dos_grid_structure * 4648* * * 4649* ******************************** 4650 subroutine control_dos_grid_structure(grid) 4651 implicit none 4652 integer grid(3) 4653 4654#include "bafdecls.fh" 4655#include "btdb.fh" 4656 4657* **** control_rtdb common block **** 4658 integer rtdb 4659 common / control_rtdb1 / rtdb 4660 4661 if (.not.btdb_get(rtdb,'band:dos-grid',mt_int,3,grid)) then 4662 grid(1) = 1 4663 grid(2) = 1 4664 grid(3) = 1 4665 end if 4666 4667 return 4668 end 4669 4670 4671 4672 4673* *********************************** 4674* * * 4675* * control_reset_band_structure * 4676* * * 4677* *********************************** 4678 subroutine control_reset_band_structure() 4679 implicit none 4680 4681#include "bafdecls.fh" 4682#include "btdb.fh" 4683#include "control.fh" 4684 4685* **** control_rtdb common block **** 4686 integer rtdb 4687 common / control_rtdb1 / rtdb 4688 4689 if (.not. btdb_cget(rtdb, 'pspw:input bvectors', 4690 > 1,input_wavefunction_filename)) 4691 > input_wavefunction_filename = 'atomic' 4692 4693 if (.not. btdb_cget(rtdb, 'pspw:output bvectors', 4694 > 1,output_wavefunction_filename)) 4695 > output_wavefunction_filename = ' ' 4696 if (output_wavefunction_filename.eq.' ')then 4697 if (input_wavefunction_filename.eq.'atomic')then 4698 call util_file_prefix('bmovecs',output_wavefunction_filename) 4699 else 4700 output_wavefunction_filename = input_wavefunction_filename 4701 endif 4702 endif 4703 if (input_wavefunction_filename.eq.'atomic')then 4704 input_wavefunction_filename = output_wavefunction_filename 4705 end if 4706 4707 4708 return 4709 end 4710 4711* ************************************** 4712* * * 4713* * control_num_kvectors_structure * 4714* * * 4715* ************************************** 4716 integer function control_num_kvectors_structure() 4717 implicit none 4718 4719#include "bafdecls.fh" 4720#include "btdb.fh" 4721#include "errquit.fh" 4722 4723* **** control_rtdb common block **** 4724 integer rtdb 4725 common / control_rtdb1 / rtdb 4726 4727* **** local variables **** 4728 logical value 4729 character*50 zone_name 4730 character*50 rtdb_name 4731 integer num_kvectors,l 4732 4733 value = btdb_cget(rtdb,'band_structure:zone_name',1,zone_name) 4734 4735 l = index(zone_name,' ') -1 4736 rtdb_name = zone_name(1:l)//':number_kvectors' 4737 value = value.and. 4738 > btdb_get(rtdb,rtdb_name,mt_int,1,num_kvectors) 4739 4740 if (.not. value) 4741 > call errquit('control_num_kvectors_structure: failed', 4742 > 0, RTDB_ERR) 4743 4744 control_num_kvectors_structure = num_kvectors 4745 return 4746 end 4747 4748 4749* ************************************ 4750* * * 4751* * control_ksvector_structure * 4752* * * 4753* ************************************ 4754 subroutine control_ksvector_structure(i,ks) 4755 implicit none 4756 integer i 4757 real*8 ks(4) 4758 4759#include "bafdecls.fh" 4760#include "btdb.fh" 4761#include "errquit.fh" 4762 4763* **** control_rtdb common block **** 4764 integer rtdb 4765 common / control_rtdb1 / rtdb 4766 4767* **** local variables **** 4768 logical value 4769 character*50 zone_name 4770 character*50 rtdb_name 4771 integer num_kvectors,l 4772 integer kvs(2) 4773 4774* **** external functions **** 4775 integer control_num_kvectors_structure 4776 external control_num_kvectors_structure 4777 4778 num_kvectors = control_num_kvectors_structure() 4779 value = BA_push_get(mt_dbl,(4*num_kvectors), 4780 > 'kvs',kvs(2),kvs(1)) 4781 if (.not. value) 4782 > call errquit('control_ksvector: failed to get zone name', 0, 4783 & MA_ERR) 4784 4785 value = value.and. 4786 > btdb_cget(rtdb,'band_structure:zone_name',1,zone_name) 4787 if (.not. value) 4788 > call errquit('control_ksvector: failed to get zone name', 0, 4789 & RTDB_ERR) 4790 4791 l = index(zone_name,' ') -1 4792 rtdb_name = zone_name(1:l)//':kvectors' 4793 value = value.and. 4794 > btdb_get(rtdb,rtdb_name,mt_dbl, 4795 > (4*num_kvectors), 4796 > dbl_mb(kvs(1))) 4797 4798 if (.not. value) 4799 > call errquit('control_ksvector: failed to get kvs', 0, 4800 & RTDB_ERR) 4801 4802 ks(1) = dbl_mb(kvs(1)+4*(i-1)) 4803 ks(2) = dbl_mb(kvs(1)+4*(i-1)+1) 4804 ks(3) = dbl_mb(kvs(1)+4*(i-1)+2) 4805 ks(4) = dbl_mb(kvs(1)+4*(i-1)+3) 4806 4807 value = value.and.BA_pop_stack(kvs(2)) 4808 4809 if (.not. value) 4810 > call errquit('control_ksvector: failed to free stack', 0, 4811 & MA_ERR) 4812 return 4813 end 4814 4815 4816* ************************************ 4817* * * 4818* * control_kvector_structure * 4819* * * 4820* ************************************ 4821 subroutine control_kvector_structure(i,kv) 4822 implicit none 4823 integer i 4824 real*8 kv(3) 4825 4826* **** local variables **** 4827 real*8 ks(4) 4828 4829* **** external functions **** 4830 real*8 lattice_unitg 4831 external lattice_unitg 4832 4833 call control_ksvector_structure(i,ks) 4834 4835 kv(1) = ks(1)*lattice_unitg(1,1) 4836 > + ks(2)*lattice_unitg(1,2) 4837 > + ks(3)*lattice_unitg(1,3) 4838 kv(2) = ks(1)*lattice_unitg(2,1) 4839 > + ks(2)*lattice_unitg(2,2) 4840 > + ks(3)*lattice_unitg(2,3) 4841 kv(3) = ks(1)*lattice_unitg(3,1) 4842 > + ks(2)*lattice_unitg(3,2) 4843 > + ks(3)*lattice_unitg(3,3) 4844 4845 return 4846 end 4847 4848 4849 4850* ***************************** 4851* * * 4852* * control_print * 4853* * * 4854* ***************************** 4855 logical function control_print(level) 4856 implicit none 4857 integer level 4858 4859 4860* **** control_print common block **** 4861 integer print_level 4862 common / control_print1 / print_level 4863 4864 4865 logical value 4866 4867 if (level.le.print_level) then 4868 value = .true. 4869 else 4870 value = .false. 4871 end if 4872 4873 control_print = value 4874 return 4875 end 4876 4877 4878* ***************************** 4879* * * 4880* * control_reduce_print * 4881* * * 4882* ***************************** 4883 subroutine control_reduce_print() 4884 implicit none 4885 4886* **** control_print common block **** 4887 integer print_level 4888 common / control_print1 / print_level 4889 4890 print_level = print_level -10 4891 return 4892 end 4893 4894* ***************************** 4895* * * 4896* * control_up_print * 4897* * * 4898* ***************************** 4899 subroutine control_up_print() 4900 implicit none 4901 4902* **** control_print common block **** 4903 integer print_level 4904 common / control_print1 / print_level 4905 4906 print_level = print_level + 10 4907 return 4908 end 4909 4910 4911 4912 4913* ************************************ 4914* * * 4915* * control_optimize_cell_strategy * 4916* * * 4917* ************************************ 4918 4919 integer function control_optimize_cell_strategy() 4920 implicit none 4921 4922#include "bafdecls.fh" 4923#include "btdb.fh" 4924 4925* **** control_rtdb common block **** 4926 integer rtdb 4927 common / control_rtdb1 / rtdb 4928 4929 integer optimize_strategy 4930 4931 if (.not.btdb_get(rtdb,'cell_optimize:optimize_strategy', 4932 > mt_int,1,optimize_strategy)) 4933 > optimize_strategy = 0 4934 4935 control_optimize_cell_strategy = optimize_strategy 4936 return 4937 end 4938 4939 4940 4941 4942* ************************************* 4943* * * 4944* * control_optimize_lattice_vectors * 4945* * * 4946* ************************************* 4947 4948 integer function control_optimize_lattice_vectors(u,v) 4949 implicit none 4950 integer u,v 4951 4952#include "bafdecls.fh" 4953#include "btdb.fh" 4954 4955* **** control_rtdb common block **** 4956 integer rtdb 4957 common / control_rtdb1 / rtdb 4958 4959 integer optimize_lattice_vectors(3,3) 4960 4961 if (.not.btdb_get(rtdb,'cell_optimize:optimize_lattice_vectors', 4962 > mt_int,9,optimize_lattice_vectors)) then 4963 4964 optimize_lattice_vectors(1,1) =1 4965 optimize_lattice_vectors(2,1) =1 4966 optimize_lattice_vectors(3,1) =1 4967 optimize_lattice_vectors(1,2) =1 4968 optimize_lattice_vectors(2,2) =1 4969 optimize_lattice_vectors(3,2) =1 4970 optimize_lattice_vectors(1,3) =1 4971 optimize_lattice_vectors(2,3) =1 4972 optimize_lattice_vectors(3,3) =1 4973 end if 4974 4975 control_optimize_lattice_vectors = optimize_lattice_vectors(u,v) 4976 return 4977 end 4978 4979* ************************************* 4980* * * 4981* * control_optimize_lattice * 4982* * * 4983* ************************************* 4984 4985 integer function control_optimize_lattice(i) 4986 implicit none 4987 integer i 4988 4989#include "bafdecls.fh" 4990#include "btdb.fh" 4991 4992* **** control_rtdb common block **** 4993 integer rtdb 4994 common / control_rtdb1 / rtdb 4995 4996 integer optimize_lattice(6) 4997 4998 if (.not.btdb_get(rtdb,'cell_optimize:optimize_lattice', 4999 > mt_int,6,optimize_lattice)) then 5000 5001 optimize_lattice(1) =1 5002 optimize_lattice(2) =1 5003 optimize_lattice(3) =1 5004 optimize_lattice(4) =1 5005 optimize_lattice(5) =1 5006 optimize_lattice(6) =1 5007 end if 5008 5009 control_optimize_lattice = optimize_lattice(i) 5010 return 5011 end 5012 5013 5014 5015* ************************************ 5016* * * 5017* * control_lmax_multipole * 5018* * * 5019* ************************************ 5020 5021 integer function control_lmax_multipole() 5022 implicit none 5023 5024#include "bafdecls.fh" 5025#include "btdb.fh" 5026 5027* **** control_rtdb common block **** 5028 integer rtdb 5029 common / control_rtdb1 / rtdb 5030 5031 integer lmax 5032 5033 if (.not.btdb_get(rtdb,'nwpw:lmax_multipole',mt_int,1,lmax)) 5034 > lmax = 0 5035 5036 control_lmax_multipole = lmax 5037 return 5038 end 5039 5040 5041 5042* ************************************ 5043* * * 5044* * control_pfft3_qsize * 5045* * * 5046* ************************************ 5047 5048 integer function control_pfft3_qsize() 5049 implicit none 5050 5051#include "bafdecls.fh" 5052#include "btdb.fh" 5053 5054* **** control_rtdb common block **** 5055 integer rtdb 5056 common / control_rtdb1 / rtdb 5057 5058 integer qmax 5059 5060 if (.not.btdb_get(rtdb,'nwpw:pfft3_qsize',mt_int,1,qmax)) 5061 > qmax = 4 5062 5063 control_pfft3_qsize = qmax 5064 return 5065 end 5066 5067 5068* ************************************ 5069* * * 5070* * control_nprj_mult * 5071* * * 5072* ************************************ 5073 5074 integer function control_nprj_mult() 5075 implicit none 5076 5077#include "bafdecls.fh" 5078#include "btdb.fh" 5079 5080* **** control_rtdb common block **** 5081 integer rtdb 5082 common / control_rtdb1 / rtdb 5083 5084 integer qmax 5085 5086 if (.not.btdb_get(rtdb,'nwpw:nprj_mult',mt_int,1,qmax)) 5087 > qmax = 1 5088 5089 control_nprj_mult = qmax 5090 return 5091 end 5092 5093 5094 5095* ************************************ 5096* * * 5097* * control_symmetry * 5098* * * 5099* ************************************ 5100 5101 integer function control_symmetry() 5102 implicit none 5103 5104#include "control.fh" 5105 5106 control_symmetry = symm_number 5107 return 5108 end 5109 5110* ************************************ 5111* * * 5112* * control_spin_orbit * 5113* * * 5114* ************************************ 5115 5116 logical function control_spin_orbit() 5117 implicit none 5118 5119#include "control.fh" 5120 5121 control_spin_orbit = spin_orbit 5122 return 5123 end 5124 5125* ************************************ 5126* * * 5127* * control_fast_erf * 5128* * * 5129* ************************************ 5130 5131 logical function control_fast_erf() 5132 implicit none 5133 5134#include "control.fh" 5135 5136 control_fast_erf = fast_erf 5137 return 5138 end 5139 5140 5141 5142* ************************************ 5143* * * 5144* * control_fmm * 5145* * * 5146* ************************************ 5147 5148 logical function control_fmm() 5149 implicit none 5150 5151#include "control.fh" 5152 5153 control_fmm = fmm 5154 return 5155 end 5156 5157* ************************************ 5158* * * 5159* * control_periodic_dipole * 5160* * * 5161* ************************************ 5162 logical function control_periodic_dipole() 5163 implicit none 5164 5165#include "control.fh" 5166 5167 control_periodic_dipole = periodic_dipole 5168 return 5169 end 5170 5171* ************************************ 5172* * * 5173* * control_smooth_cutoff * 5174* * * 5175* ************************************ 5176 logical function control_smooth_cutoff() 5177 implicit none 5178 5179#include "control.fh" 5180 5181 control_smooth_cutoff = smooth_cutoff 5182 return 5183 end 5184 5185* ************************************ 5186* * * 5187* * control_smooth_cutoff_values * 5188* * * 5189* ************************************ 5190 real*8 function control_smooth_cutoff_values(i) 5191 implicit none 5192 integer i 5193 5194#include "control.fh" 5195 5196 control_smooth_cutoff_values = smooth_cutoff_values(i) 5197 return 5198 end 5199 5200 5201 5202 5203* ************************************ 5204* * * 5205* * control_fmm_lmax * 5206* * * 5207* ************************************ 5208 5209 integer function control_fmm_lmax() 5210 implicit none 5211 5212#include "control.fh" 5213 5214 control_fmm_lmax = fmm_lmax 5215 return 5216 end 5217 5218* ************************************ 5219* * * 5220* * control_fmm_lr * 5221* * * 5222* ************************************ 5223 5224 integer function control_fmm_lr() 5225 implicit none 5226 5227#include "control.fh" 5228 5229 control_fmm_lr = fmm_lr 5230 return 5231 end 5232 5233 5234 5235 5236 5237* ***************************** 5238* * * 5239* * control_pressure * 5240* * * 5241* ***************************** 5242 logical function control_pressure() 5243 implicit none 5244 5245#include "bafdecls.fh" 5246#include "btdb.fh" 5247 5248* **** control_rtdb common block **** 5249 integer rtdb 5250 common / control_rtdb1 / rtdb 5251 5252 logical value 5253 5254 value = .false. 5255 if (.not.btdb_get(rtdb,'cpmd:pressure',mt_log,1,value)) 5256 > value = .false. 5257 5258 control_pressure = value 5259 return 5260 end 5261 5262 5263* *********************************** 5264* * * 5265* * control_init_velocities * 5266* * * 5267* *********************************** 5268 logical function control_init_velocities() 5269 implicit none 5270 5271#include "bafdecls.fh" 5272#include "btdb.fh" 5273 5274* **** control_rtdb common block **** 5275 integer rtdb 5276 common / control_rtdb1 / rtdb 5277 5278 logical ok,tmp 5279 5280 if (.not.btdb_get(rtdb,'nwpw:init_velocities',mt_log,1,tmp)) 5281 > tmp = .false. 5282 5283 if (tmp) ok = rtdb_delete(rtdb,'nwpw:init_velocities') 5284 5285 control_init_velocities = tmp 5286 return 5287 end 5288 5289 5290* *********************************** 5291* * * 5292* * control_kbpp_ray * 5293* * * 5294* *********************************** 5295 logical function control_kbpp_ray() 5296 implicit none 5297 5298#include "bafdecls.fh" 5299#include "btdb.fh" 5300 5301* **** control_rtdb common block **** 5302 integer rtdb 5303 common / control_rtdb1 / rtdb 5304 5305 logical value 5306 5307 value = .true. 5308 if (.not.btdb_get(rtdb,'nwpw:kbpp_ray', 5309 > mt_log,1,value)) 5310 > value = .false. 5311 5312 control_kbpp_ray = value 5313 return 5314 end 5315 5316 5317 5318* *********************************** 5319* * * 5320* * control_kbpp_filter * 5321* * * 5322* *********************************** 5323 logical function control_kbpp_filter() 5324 implicit none 5325 5326#include "bafdecls.fh" 5327#include "btdb.fh" 5328 5329* **** control_rtdb common block **** 5330 integer rtdb 5331 common / control_rtdb1 / rtdb 5332 5333 logical value 5334 5335 value = .true. 5336 if (.not.btdb_get(rtdb,'nwpw:kbpp_filter', 5337 > mt_log,1,value)) 5338 > value = .false. 5339 5340 control_kbpp_filter = value 5341 return 5342 end 5343 5344 5345* *********************************** 5346* * * 5347* * control_brillioun_ondisk * 5348* * * 5349* *********************************** 5350 logical function control_brillioun_ondisk() 5351 implicit none 5352 5353#include "bafdecls.fh" 5354#include "btdb.fh" 5355 5356* **** control_rtdb common block **** 5357 integer rtdb 5358 common / control_rtdb1 / rtdb 5359 5360 logical value 5361 5362 value = .true. 5363 if (.not.btdb_get(rtdb,'nwpw:brillioun_ondisk', 5364 > mt_log,1,value)) 5365 > value = .false. 5366 5367 control_brillioun_ondisk = value 5368 return 5369 end 5370 5371 5372 5373* ************************************** 5374* * * 5375* * control_brillioun_use_symmetry * 5376* * * 5377* ************************************** 5378 logical function control_brillioun_use_symmetry() 5379 implicit none 5380 5381#include "bafdecls.fh" 5382#include "btdb.fh" 5383 5384* **** control_rtdb common block **** 5385 integer rtdb 5386 common / control_rtdb1 / rtdb 5387 5388 logical value 5389 5390 value = .true. 5391 if (.not.btdb_get(rtdb,'nwpw:brillioun_use_symmetry', 5392 > mt_log,1,value)) 5393 > value = .true. 5394 5395 control_brillioun_use_symmetry = value 5396 return 5397 end 5398 5399 5400 5401 5402* *********************************** 5403* * * 5404* * control_mparallelized * 5405* * * 5406* *********************************** 5407 logical function control_mparallelized() 5408 implicit none 5409 5410#include "bafdecls.fh" 5411#include "btdb.fh" 5412 5413* **** control_rtdb common block **** 5414 integer rtdb 5415 common / control_rtdb1 / rtdb 5416 5417 logical value 5418 5419 value = .true. 5420 if (.not.btdb_get(rtdb,'nwpw:mparallelized',mt_log,1,value)) 5421 > value = .false. 5422 5423 control_mparallelized = value 5424 return 5425 end 5426 5427 5428* *********************************** 5429* * * 5430* * control_mreplicate_size * 5431* * * 5432* *********************************** 5433 integer function control_mreplicate_size() 5434 implicit none 5435 5436#include "bafdecls.fh" 5437#include "btdb.fh" 5438 5439* **** control_rtdb common block **** 5440 integer rtdb 5441 common / control_rtdb1 / rtdb 5442 5443 integer value 5444 5445 value = 16 5446 if (.not.btdb_get(rtdb,'nwpw:mreplicate_size',mt_int,1,value)) 5447 > value = 16 5448 5449 control_mreplicate_size = value 5450 return 5451 end 5452 5453 5454* *********************************** 5455* * * 5456* * control_bo_cpmd * 5457* * * 5458* *********************************** 5459 logical function control_bo_cpmd() 5460 implicit none 5461 5462#include "bafdecls.fh" 5463#include "btdb.fh" 5464 5465* **** control_rtdb common block **** 5466 integer rtdb 5467 common / control_rtdb1 / rtdb 5468 5469 logical value 5470 5471 value = .true. 5472 if (.not.btdb_get(rtdb,'nwpw:bo_cpmd',mt_log,1,value)) 5473 > value = .false. 5474 5475 control_bo_cpmd = value 5476 return 5477 end 5478 5479* *********************************** 5480* * * 5481* * control_hfxon_virtual * 5482* * * 5483* *********************************** 5484 logical function control_hfxon_virtual() 5485 implicit none 5486 5487#include "bafdecls.fh" 5488#include "btdb.fh" 5489 5490* **** control_rtdb common block **** 5491 integer rtdb 5492 common / control_rtdb1 / rtdb 5493 5494 logical value 5495 5496 value = .true. 5497 if (.not.btdb_get(rtdb,'nwpw:hfxon_virtual',mt_log,1,value)) 5498 > value = .true. 5499 5500 control_hfxon_virtual = value 5501 return 5502 end 5503 5504 5505* *********************************** 5506* * * 5507* * control_gradient_virtual * 5508* * * 5509* *********************************** 5510 logical function control_gradient_virtual() 5511 implicit none 5512 5513#include "bafdecls.fh" 5514#include "btdb.fh" 5515 5516* **** control_rtdb common block **** 5517 integer rtdb 5518 common / control_rtdb1 / rtdb 5519 5520 logical value 5521 5522 value = .false. 5523 if (.not.btdb_get(rtdb,'nwpw:gradient_virtual',mt_log,1,value)) 5524 > value = .false. 5525 5526 control_gradient_virtual = value 5527 return 5528 end 5529 5530* *********************************** 5531* * * 5532* * control_gradient_virtual_fac * 5533* * * 5534* *********************************** 5535 integer function control_gradient_virtual_fac() 5536 implicit none 5537 5538#include "bafdecls.fh" 5539#include "btdb.fh" 5540 5541* **** control_rtdb common block **** 5542 integer rtdb 5543 common / control_rtdb1 / rtdb 5544 5545 integer fac 5546 5547 fac = 19 5548 if (.not.btdb_get(rtdb,'nwpw:gradient_virtual_fac',mt_log,1,fac)) 5549 > fac = 19 5550 5551 control_gradient_virtual_fac = fac 5552 return 5553 end 5554 5555 5556* *********************************** 5557* * * 5558* * control_bound_virtual * 5559* * * 5560* *********************************** 5561 logical function control_bound_virtual() 5562 implicit none 5563 5564#include "bafdecls.fh" 5565#include "btdb.fh" 5566 5567* **** control_rtdb common block **** 5568 integer rtdb 5569 common / control_rtdb1 / rtdb 5570 5571 logical value 5572 5573 value = .false. 5574 if (.not.btdb_get(rtdb,'nwpw:bound_virtual',mt_log,1,value)) 5575 > value = .false. 5576 5577 control_bound_virtual = value 5578 return 5579 end 5580 5581* *********************************** 5582* * * 5583* * control_confine_virtual * 5584* * * 5585* *********************************** 5586 logical function control_confine_virtual() 5587 implicit none 5588 5589#include "bafdecls.fh" 5590#include "btdb.fh" 5591 5592* **** control_rtdb common block **** 5593 integer rtdb 5594 common / control_rtdb1 / rtdb 5595 5596 logical value 5597 5598 value = .false. 5599 if (.not.btdb_get(rtdb,'nwpw:confine_virtual',mt_log,1,value)) 5600 > value = .false. 5601 5602 control_confine_virtual = value 5603 return 5604 end 5605 5606 5607 5608 5609 5610 5611 5612 5613* ********************************************** 5614* * * 5615* * control_init_velocities_temperature * 5616* * * 5617* ********************************************** 5618 real*8 function control_init_velocities_temperature() 5619 implicit none 5620 5621#include "bafdecls.fh" 5622#include "btdb.fh" 5623 5624* **** control_rtdb common block **** 5625 integer rtdb 5626 common / control_rtdb1 / rtdb 5627 5628 real*8 value 5629 5630 if (.not.btdb_get(rtdb,'nwpw:init_velocities_temperature', 5631 > mt_dbl,1,value)) 5632 > value = 300.0d0 5633 5634 control_init_velocities_temperature = value 5635 return 5636 end 5637 5638* ********************************************** 5639* * * 5640* * control_init_velocities_seed * 5641* * * 5642* ********************************************** 5643 integer function control_init_velocities_seed() 5644 implicit none 5645 5646#include "bafdecls.fh" 5647#include "btdb.fh" 5648 5649* **** control_rtdb common block **** 5650 integer rtdb 5651 common / control_rtdb1 / rtdb 5652 5653 integer seed 5654 5655 if (.not.btdb_get(rtdb,'nwpw:init_velocities_seed', 5656 > mt_int,1,seed)) 5657 > seed = 494 5658 5659 control_init_velocities_seed = seed 5660 return 5661 end 5662 5663 5664 5665 5666* *********************************** 5667* * * 5668* * control_wannier_timestep * 5669* * * 5670* *********************************** 5671 real*8 function control_wannier_timestep() 5672 implicit none 5673 5674#include "bafdecls.fh" 5675#include "btdb.fh" 5676 5677* **** control_rtdb common block **** 5678 integer rtdb 5679 common / control_rtdb1 / rtdb 5680 5681 real*8 value 5682 5683 if (.not.btdb_get(rtdb,'wannier:time_step',mt_dbl,1,value)) 5684 > value = 2.7e-2 5685 5686 control_wannier_timestep = value 5687 return 5688 end 5689 5690 5691 5692* *********************************** 5693* * * 5694* * control_wannier_maxiter * 5695* * * 5696* *********************************** 5697 integer function control_wannier_maxiter() 5698 implicit none 5699 5700#include "bafdecls.fh" 5701#include "btdb.fh" 5702 5703* **** control_rtdb common block **** 5704 integer rtdb 5705 common / control_rtdb1 / rtdb 5706 5707 integer value 5708 5709 if (.not.btdb_get(rtdb,'wannier:maxiter',mt_int,1,value)) 5710 > value = 500 5711 5712 control_wannier_maxiter = value 5713 return 5714 end 5715 5716 5717* ********************************************** 5718* * * 5719* * control_attenuation * 5720* * * 5721* ********************************************** 5722 real*8 function control_attenuation() 5723 implicit none 5724#include "control.fh" 5725 control_attenuation = attenuation 5726 return 5727 end 5728 5729 5730 5731ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 5732 subroutine set_two_component_pseudopotential() 5733ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 5734c this called to signal that a two_component ppot is in use 5735ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 5736#include "control.fh" 5737 two_comp_ppot=.true. 5738 return 5739 end 5740ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 5741 logical function two_component_pseudopotential() 5742#include "control.fh" 5743 two_component_pseudopotential=two_comp_ppot 5744 return 5745 end 5746ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 5747 5748 5749 5750 5751* *********************************** 5752* * * 5753* * control_ecut_wcut_default * 5754* * * 5755* *********************************** 5756 subroutine control_ecut_wcut_default(rtdb,ecut,wcut) 5757 implicit none 5758 integer rtdb 5759 real*8 ecut,wcut 5760 5761#include "bafdecls.fh" 5762#include "btdb.fh" 5763#include "beom.fh" 5764#include "errquit.fh" 5765 5766* **** local variables **** 5767 integer taskid,MASTER 5768 parameter (MASTER=0) 5769 5770 logical value 5771 integer geom,ii,nion,l,h 5772 real*8 wcut_max,wcut_atom,q,rion(3) 5773 character*4 aname 5774 character*16 t,rtdbname 5775 character*255 sdir_name 5776 5777 5778* **** external functions **** 5779 real*8 nwpw_libgcutoff 5780 external nwpw_libgcutoff 5781 character*4 ion_aname_geom 5782 external ion_aname_geom 5783 logical parsepointcharge 5784 external parsepointcharge 5785 5786 !*** parse psp information on rtdb *** 5787 if (.not.btdb_get(rtdb,'nwpw:psp:cutoff',mt_dbl,1,wcut_max)) then 5788 call Parallel_taskid(taskid) 5789 wcut_max = -1.0d0 5790 value = beom_create(geom,'geometry') 5791 value = value.and.beom_rtdb_load(rtdb,geom,'geometry') 5792 value = value.and.geom_ncent(geom,nion) 5793 if (.not. value) 5794 > call errquit('control_ecut_wcut_default:cannot load geometry', 5795 > 0,GEOM_ERR) 5796 5797 do ii=1,nion 5798 if (.not.geom_cent_get(geom,ii,t,rion,q)) 5799 > call errquit('control_ecut_wcut_default:error reading ions', 5800 > 0,GEOM_ERR) 5801 if (.not.parsepointcharge(t)) then 5802 aname = ion_aname_geom(geom,ii) 5803 l = index(aname,' ') - 1 5804 if (l.le.0) l = 4 5805 rtdbname = aname(1:l)//':cutoff' 5806 l = l+7 5807 if (.not.btdb_get(rtdb,rtdbname(1:l),mt_dbl,1,wcut_atom)) 5808 > then 5809 5810 !*** define wcut_atom from psplibrary *** 5811 value = btdb_parallel(.false.) 5812 if (taskid.eq.MASTER) then 5813 call util_directory_name(sdir_name,.true.,0) 5814 h = index(sdir_name,' ') - 1 5815 open(unit=99,file=sdir_name(1:h)//'/junk.inp', 5816 > status='unknown') 5817 close(unit=99,status='delete') 5818 5819 call nwpw_libgeninp(1,aname, 5820 > sdir_name(1:h)//'/junk.inp') 5821 wcut_atom=nwpw_libgcutoff(aname) 5822 5823 if(.not.btdb_put(rtdb,rtdbname(1:l), 5824 > mt_dbl,1,wcut_atom)) 5825 > call errquit( 5826 > 'control_ecut_wcut_default:cannot write wcut_atom',0,0) 5827 5828 end if 5829 value = btdb_parallel(.true.) 5830 call Parallel_Brdcst_value(MASTER,wcut_atom) 5831 5832 end if 5833 5834 !*** reset wcut_max *** 5835 if (wcut_atom.gt.wcut_max) then 5836 wcut_max = wcut_atom 5837 value = btdb_parallel(.false.) 5838 if (taskid.eq.MASTER) then 5839 if (.not.btdb_put(rtdb,'nwpw:psp:cutoff', 5840 > mt_dbl,1,wcut_max)) 5841 > call errquit( 5842 > 'control_ecut_wcut_default:cannot write wcut_max',0,0) 5843 end if 5844 value = btdb_parallel(.true.) 5845 end if 5846 end if 5847 end do 5848 if (.not. beom_destroy(geom)) 5849 > call errquit('control_ecut_wcut_default:cannot destroy geom', 5850 > 0,GEOM_ERR) 5851 end if 5852 5853 wcut = wcut_max 5854 ecut = 2.0d0*wcut 5855 5856 return 5857 end 5858 5859 5860 5861* *********************************** 5862* * * 5863* * control_pspspin * 5864* * * 5865* *********************************** 5866 logical function control_pspspin() 5867 implicit none 5868 5869#include "bafdecls.fh" 5870#include "btdb.fh" 5871 5872* **** control_rtdb common block **** 5873 integer rtdb 5874 common / control_rtdb1 / rtdb 5875 5876 logical value 5877 5878 value = .true. 5879 if (.not.btdb_get(rtdb,'nwpw:pspspin',mt_log,1,value)) 5880 > value = .false. 5881 5882 control_pspspin = value 5883 return 5884 end 5885 5886 5887* *********************************** 5888* * * 5889* * control_set_pspspin * 5890* * * 5891* *********************************** 5892 subroutine control_set_pspspin(nion, 5893 > upscale,downscale, 5894 > upl,downl, 5895 > upm,downm, 5896 > upions,downions) 5897 implicit none 5898 integer nion 5899 real*8 upscale(*),downscale(*) 5900 integer upl(*),downl(*),upm(*),downm(*) 5901 logical upions(*),downions(*) 5902 5903#include "bafdecls.fh" 5904#include "btdb.fh" 5905#include "errquit.fh" 5906 5907* **** control_rtdb common block **** 5908 integer rtdb 5909 common / control_rtdb1 / rtdb 5910 5911 integer taskid,MASTER 5912 parameter (MASTER=0) 5913 logical iamup 5914 integer i,ma_type,nactive_atoms,h_actlist,l_actlist 5915 integer pcount,ip,tl,tm 5916 real*8 tscale 5917 character*50 rtdb_name 5918 5919* *** external functions *** 5920 character*7 c_index_name 5921 external c_index_name 5922 5923 call Parallel_taskid(taskid) 5924 if (.not.btdb_get(rtdb,'nwpw:pspspin_count',mt_int,1,pcount)) 5925 > pcount = 0 5926 5927 if (taskid.eq.MASTER) write(*,2300) 5928 do ip=1,pcount 5929 5930 rtdb_name = 'nwpw:pspspin_iamup:'//c_index_name(ip) 5931 if (.not.btdb_get(rtdb,rtdb_name,mt_log,1,iamup)) iamup=.true. 5932 5933 if (iamup) then 5934 rtdb_name = 'nwpw:pspspin_upscale:'//c_index_name(ip) 5935 if (.not.btdb_get(rtdb,rtdb_name,mt_dbl,1,tscale)) 5936 > tscale = 1.0d0 5937 rtdb_name = 'nwpw:pspspin_upl:'//c_index_name(ip) 5938 if (.not.btdb_get(rtdb,rtdb_name,mt_int,1,tl)) 5939 > tl = -1 5940 rtdb_name = 'nwpw:pspspin_upm:'//c_index_name(ip) 5941 if (.not.btdb_get(rtdb,rtdb_name,mt_int,1,tm)) 5942 > tm = 99999 5943 if (taskid.eq.MASTER) then 5944 if (tm.lt.999) then 5945 write(*,2305) tl,tm,tscale 5946 else 5947 write(*,2301) tl,tscale 5948 end if 5949 end if 5950 5951 rtdb_name = 'nwpw:pspspin_upions:'//c_index_name(ip) 5952 if (rtdb_ma_get(rtdb,rtdb_name, ma_type, 5953 > nactive_atoms, h_actlist)) then 5954 if (.not.BA_get_index(h_actlist,l_actlist)) 5955 > call errquit( 5956 > 'control_set_pspspin: ma_get_index failed',0, 5957 > MA_ERR) 5958 5959 if (taskid.eq.MASTER) 5960 > write(*,2302) (int_mb(l_actlist+i),i=0,nactive_atoms-1) 5961 5962 do i=1,nactive_atoms 5963 upions(int_mb(l_actlist+i-1)) = .true. 5964 upl(int_mb(l_actlist+i-1)) = tl 5965 upm(int_mb(l_actlist+i-1)) = tm 5966 upscale(int_mb(l_actlist+i-1)) = tscale 5967 end do 5968 if (.not. BA_free_heap(h_actlist)) 5969 > call errquit( 5970 > 'control_set_pspspin:error freeing heap memory',0,MA_ERR) 5971 end if 5972 5973 5974 else 5975 rtdb_name = 'nwpw:pspspin_downscale:'//c_index_name(ip) 5976 if (.not.btdb_get(rtdb,rtdb_name,mt_dbl,1,tscale)) 5977 > tscale = 1.0d0 5978 rtdb_name = 'nwpw:pspspin_downl:'//c_index_name(ip) 5979 if (.not.btdb_get(rtdb,rtdb_name,mt_int,1,tl)) 5980 > tl = -1 5981 rtdb_name = 'nwpw:pspspin_downm:'//c_index_name(ip) 5982 if (.not.btdb_get(rtdb,rtdb_name,mt_int,1,tm)) 5983 > tm = 99999 5984 if (taskid.eq.MASTER) then 5985 if (tm.lt.999) then 5986 write(*,2306) tl,tm,tscale 5987 else 5988 write(*,2303) tl,tscale 5989 end if 5990 end if 5991 5992 rtdb_name = 'nwpw:pspspin_downions:'//c_index_name(ip) 5993 if (rtdb_ma_get(rtdb,rtdb_name,ma_type, 5994 > nactive_atoms,h_actlist)) then 5995 if (.not.BA_get_index(h_actlist,l_actlist)) 5996 > call errquit( 5997 > 'control_set_pspspin: ma_get_index failed',1, 5998 > MA_ERR) 5999 6000 if (taskid.eq.MASTER) 6001 > write(*,2304) (int_mb(l_actlist+i),i=0,nactive_atoms-1) 6002 6003 do i=1,nactive_atoms 6004 downions(int_mb(l_actlist+i-1)) = .true. 6005 downl(int_mb(l_actlist+i-1)) = tl 6006 downm(int_mb(l_actlist+i-1)) = tm 6007 downscale(int_mb(l_actlist+i-1)) = tscale 6008 end do 6009 if (.not.BA_free_heap(h_actlist)) 6010 > call errquit( 6011 > 'control_set_pspspin:error freeing heap memory',1,MA_ERR) 6012 end if 6013 end if 6014 end do 6015 6016 2300 format(/1x,"Antiferromagnetic Pentalty Function Input:") 6017 2301 format(2x," - pspspin: up l =",I2,10x," scale =",F8.3) 6018 2302 format(2x," - pspspin: up ion indexes =",10I5) 6019 2303 format(2x," - pspspin: down l =",I2,10x," scale =",F8.3) 6020 2304 format(2x," - pspspin: down ion indexes =",10I5) 6021 2305 format(2x," - pspspin: up l =",I2," not_m=",I3," scale =",F8.3) 6022 2306 format(2x," - pspspin: down l =",I2," not_m=",I3," scale =",F8.3) 6023 6024 return 6025 end 6026 6027* *********************************** 6028* * * 6029* * control_psputerm * 6030* * * 6031* *********************************** 6032 logical function control_psputerm() 6033 implicit none 6034 6035#include "bafdecls.fh" 6036#include "btdb.fh" 6037 6038* **** control_rtdb common block **** 6039 integer rtdb 6040 common / control_rtdb1 / rtdb 6041 6042 logical value 6043 6044 value = .true. 6045 if (.not.btdb_get(rtdb,'nwpw:uterm',mt_log,1,value)) 6046 > value = .false. 6047 6048 control_psputerm = value 6049 return 6050 end 6051 6052* *********************************** 6053* * * 6054* * control_pspnuterms * 6055* * * 6056* *********************************** 6057 integer function control_pspnuterms() 6058 implicit none 6059 6060#include "bafdecls.fh" 6061#include "btdb.fh" 6062 6063* **** control_rtdb common block **** 6064 integer rtdb 6065 common / control_rtdb1 / rtdb 6066 6067 integer nuterms 6068 6069 nuterms = 0 6070 if (.not.btdb_get(rtdb,'nwpw:nuterms',mt_int,1,nuterms)) 6071 > nuterms = 0 6072 6073 control_pspnuterms = nuterms 6074 return 6075 end 6076 6077 6078* *********************************** 6079* * * 6080* * control_set_psputerm * 6081* * * 6082* *********************************** 6083 subroutine control_set_psputerm(nion,nuterms,l,uscale,jscale,ions) 6084 implicit none 6085 integer nion,nuterms,l(nuterms) 6086 real*8 uscale(nuterms),jscale(nuterms) 6087 logical ions(nion,nuterms) 6088 6089#include "bafdecls.fh" 6090#include "btdb.fh" 6091#include "errquit.fh" 6092 6093* **** control_rtdb common block **** 6094 integer rtdb 6095 common / control_rtdb1 / rtdb 6096 6097 integer taskid,MASTER 6098 parameter (MASTER=0) 6099 integer nu,i,ma_type,nactive_atoms,h_actlist,l_actlist 6100 character*50 rtdb_name 6101 6102* **** external functions **** 6103 character*7 c_index_name 6104 external c_index_name 6105 6106 call Parallel_taskid(taskid) 6107 6108 if (taskid.eq.MASTER) write(*,3300) 6109 6110 do nu=1,nuterms 6111 6112 rtdb_name = 'nwpw:uterm_scale:'//c_index_name(nu) 6113 if (.not.btdb_get(rtdb,rtdb_name,mt_dbl,1,uscale(nu))) 6114 > uscale(nu) = 0.0d0 6115 6116 rtdb_name = 'nwpw:jterm_scale:'//c_index_name(nu) 6117 if (.not.btdb_get(rtdb,rtdb_name,mt_dbl,1,jscale(nu))) 6118 > jscale(nu) = 0.0d0 6119 6120 rtdb_name = 'nwpw:uterm_l:'//c_index_name(nu) 6121 if (.not.btdb_get(rtdb,rtdb_name,mt_int,1,l(nu))) 6122 > l(nu) = -1 6123 6124 do i=1,nion 6125 ions(i,nu) = .false. 6126 end do 6127 6128 if (taskid.eq.MASTER) write(*,3301) l(nu),uscale(nu),jscale(nu) 6129 6130 rtdb_name = 'nwpw:uterm_ions:'//c_index_name(nu) 6131 if (rtdb_ma_get(rtdb, rtdb_name,ma_type, 6132 > nactive_atoms, h_actlist)) then 6133 if (.not.BA_get_index(h_actlist,l_actlist)) 6134 > call errquit( 6135 > 'control_set_psputerm: ma_get_index failed',0, 6136 > MA_ERR) 6137 6138 if (taskid.eq.MASTER) then 6139 write(*,3302) (int_mb(l_actlist+i),i=0,nactive_atoms-1) 6140 write(*,*) 6141 end if 6142 6143 do i=1,nactive_atoms 6144 ions(int_mb(l_actlist+i-1),nu) = .true. 6145 end do 6146 if (.not. BA_free_heap(h_actlist)) 6147 > call errquit( 6148 > 'control_set_psputerm:error freeing heap memory',0,MA_ERR) 6149 end if 6150 6151 end do 6152 6153 3300 format(/1x,"Hubbard Uterm Function Input:") 6154 3301 format(2x," - uterm: l =",I2," U=",F8.3," J=",F8.3) 6155 3302 format(2x," - uterm: ion indexes =",10I5) 6156 return 6157 end 6158 6159 6160 6161* *********************************** 6162* * * 6163* * control_use_grid_cmp * 6164* * * 6165* *********************************** 6166 logical function control_use_grid_cmp() 6167 implicit none 6168 6169#include "bafdecls.fh" 6170#include "btdb.fh" 6171 6172* **** control_rtdb common block **** 6173 integer rtdb 6174 common / control_rtdb1 / rtdb 6175 6176 logical value 6177 6178 value = .true. 6179 if (.not.btdb_get(rtdb,'nwpw:use_grid_cmp',mt_log,1,value)) 6180 > value = .false. 6181 6182 control_use_grid_cmp = value 6183 return 6184 end 6185 6186 6187* *********************************** 6188* * * 6189* * control_wgc_alphabetalambda * 6190* * * 6191* *********************************** 6192 real*8 function control_wgc_alphabetalambda(i) 6193 implicit none 6194 integer i 6195 6196#include "bafdecls.fh" 6197#include "btdb.fh" 6198 6199* **** control_rtdb common block **** 6200 integer rtdb 6201 common / control_rtdb1 / rtdb 6202 6203 real*8 value 6204 6205 if (i.eq.1) then 6206 if (.not.btdb_get(rtdb,'nwpw:wgc_alpha',mt_dbl,1,value)) 6207 > value = 5.0d0/6.0d0 6208 else if (i.eq.2) then 6209 if (.not.btdb_get(rtdb,'nwpw:wgc_beta',mt_dbl,1,value)) 6210 > value = 5.0d0/6.0d0 6211 else 6212 if (.not.btdb_get(rtdb,'nwpw:wgc_lambda',mt_dbl,1,value)) 6213 > value = 1.0d0 6214 end if 6215 6216 control_wgc_alphabetalambda = value 6217 return 6218 end 6219 6220 6221 6222* *********************************** 6223* * * 6224* * control_mc_step_size * 6225* * * 6226* *********************************** 6227 real*8 function control_mc_step_size() 6228 implicit none 6229 6230#include "bafdecls.fh" 6231#include "btdb.fh" 6232 6233* **** control_rtdb common block **** 6234 integer rtdb 6235 common / control_rtdb1 / rtdb 6236 6237 real*8 value 6238 6239 if (.not.btdb_get(rtdb,'nwpw:mc_step_size',mt_dbl,1,value)) 6240 > value = 0.50d0 6241 6242 control_mc_step_size = value 6243 return 6244 end 6245 6246 6247 6248 6249* *********************************** 6250* * * 6251* * control_mc_volume_step * 6252* * * 6253* *********************************** 6254 real*8 function control_mc_volume_step() 6255 implicit none 6256 6257#include "bafdecls.fh" 6258#include "btdb.fh" 6259 6260* **** control_rtdb common block **** 6261 integer rtdb 6262 common / control_rtdb1 / rtdb 6263 6264 real*8 value 6265 6266 if (.not.btdb_get(rtdb,'nwpw:mc_volume_step',mt_dbl,1,value)) 6267 > value = 0.10d0 6268 6269 control_mc_volume_step = value 6270 return 6271 end 6272 6273 6274 6275 6276* *********************************** 6277* * * 6278* * control_mc_aratio * 6279* * * 6280* *********************************** 6281 real*8 function control_mc_aratio() 6282 implicit none 6283 6284#include "bafdecls.fh" 6285#include "btdb.fh" 6286 6287* **** control_rtdb common block **** 6288 integer rtdb 6289 common / control_rtdb1 / rtdb 6290 6291 real*8 value 6292 6293 if (.not.btdb_get(rtdb,'nwpw:mc_aratio',mt_dbl,1,value)) 6294 > value = 0.234d0 6295 6296 control_mc_aratio = value 6297 return 6298 end 6299 6300 6301* *********************************** 6302* * * 6303* * control_mc_Temperature * 6304* * * 6305* *********************************** 6306 real*8 function control_mc_Temperature() 6307 implicit none 6308 6309#include "bafdecls.fh" 6310#include "btdb.fh" 6311 6312* **** control_rtdb common block **** 6313 integer rtdb 6314 common / control_rtdb1 / rtdb 6315 6316 real*8 value 6317 6318 if (.not.btdb_get(rtdb,'nwpw:mc_temperature',mt_dbl,1,value)) 6319 > value = 298.15d0 6320 6321 control_mc_Temperature = value 6322 return 6323 end 6324 6325 6326* *********************************** 6327* * * 6328* * control_mc_pressure * 6329* * * 6330* *********************************** 6331 real*8 function control_mc_pressure() 6332 implicit none 6333 6334#include "bafdecls.fh" 6335#include "btdb.fh" 6336 6337* **** local variables **** 6338 real*8 autoMbar,autoGPa,autoatm 6339 parameter (autoMbar=294.214239071d0) 6340 parameter (autoGPa=autoMbar*100.0d0) 6341 parameter (autoatm =290.360032539d6) 6342 6343 6344* **** control_rtdb common block **** 6345 integer rtdb 6346 common / control_rtdb1 / rtdb 6347 6348 real*8 value 6349 6350 if (.not.btdb_get(rtdb,'nwpw:mc_pressure',mt_dbl,1,value)) 6351 > value = 1.0d0/autoatm 6352 6353 control_mc_pressure = value 6354 return 6355 end 6356 6357 6358* *********************************** 6359* * * 6360* * control_mc_ddx * 6361* * * 6362* *********************************** 6363 real*8 function control_mc_ddx() 6364 implicit none 6365 6366#include "bafdecls.fh" 6367#include "btdb.fh" 6368 6369* **** control_rtdb common block **** 6370 integer rtdb 6371 common / control_rtdb1 / rtdb 6372 6373 real*8 value 6374 6375 if (.not.btdb_get(rtdb,'nwpw:mc_ddx',mt_dbl,1,value)) 6376 > value = 0.0d0 6377 6378 control_mc_ddx = value 6379 return 6380 end 6381 6382 6383* *********************************** 6384* * * 6385* * control_mc_ddv * 6386* * * 6387* *********************************** 6388 real*8 function control_mc_ddv() 6389 implicit none 6390 6391#include "bafdecls.fh" 6392#include "btdb.fh" 6393 6394* **** control_rtdb common block **** 6395 integer rtdb 6396 common / control_rtdb1 / rtdb 6397 6398 real*8 value 6399 6400 if (.not.btdb_get(rtdb,'nwpw:mc_ddv',mt_dbl,1,value)) 6401 > value = 0.0d0 6402 6403 control_mc_ddv = value 6404 return 6405 end 6406 6407 6408 6409* *********************************** 6410* * * 6411* * control_mc_seed * 6412* * * 6413* *********************************** 6414 integer function control_mc_seed() 6415 implicit none 6416 6417#include "bafdecls.fh" 6418#include "btdb.fh" 6419 6420* **** control_rtdb common block **** 6421 integer rtdb 6422 common / control_rtdb1 / rtdb 6423 6424 integer ivalue 6425 6426 if (.not.btdb_get(rtdb,'nwpw:mc_seed',mt_int,1,ivalue)) 6427 > ivalue = 9484943 6428 6429 control_mc_seed = ivalue 6430 return 6431 end 6432 6433 6434* *********************************** 6435* * * 6436* * control_mc_algorithm * 6437* * * 6438* *********************************** 6439 integer function control_mc_algorithm() 6440 implicit none 6441 6442#include "bafdecls.fh" 6443#include "btdb.fh" 6444 6445* **** control_rtdb common block **** 6446 integer rtdb 6447 common / control_rtdb1 / rtdb 6448 6449 integer ivalue 6450 6451 if (.not.btdb_get(rtdb,'nwpw:mc_algorithm',mt_int,1,ivalue)) 6452 > ivalue = 1 6453 6454 control_mc_algorithm = ivalue 6455 return 6456 end 6457 6458 6459 6460 6461* *********************************** 6462* * * 6463* * control_mc_atom_direction * 6464* * * 6465* *********************************** 6466 subroutine control_mc_atom_direction(mc_atom_direction) 6467 implicit none 6468 real*8 mc_atom_direction(3) 6469 6470#include "bafdecls.fh" 6471#include "btdb.fh" 6472 6473* **** control_rtdb common block **** 6474 integer rtdb 6475 common / control_rtdb1 / rtdb 6476 6477 if (.not.btdb_get(rtdb,'nwpw:mc_atom_direction', 6478 > mt_dbl,3,mc_atom_direction)) then 6479 mc_atom_direction(1) = 1.0d0 6480 mc_atom_direction(2) = 1.0d0 6481 mc_atom_direction(3) = 1.0d0 6482 end if 6483 6484 return 6485 end 6486 6487 6488 6489* *********************************** 6490* * * 6491* * control_mc_ngroups * 6492* * * 6493* *********************************** 6494 subroutine control_mc_ngroups(mc_napply,mc_ngroups,mc_group_size) 6495 implicit none 6496 integer mc_napply,mc_ngroups,mc_group_size 6497 6498#include "bafdecls.fh" 6499#include "btdb.fh" 6500 6501* **** control_rtdb common block **** 6502 integer rtdb 6503 common / control_rtdb1 / rtdb 6504 6505 if (.not.btdb_get(rtdb,'nwpw:mc_napply', 6506 > mt_int,1,mc_napply)) then 6507 mc_napply = 1 6508 end if 6509 if (.not.btdb_get(rtdb,'nwpw:mc_ngroups', 6510 > mt_int,1,mc_ngroups)) then 6511 mc_ngroups = 0 6512 end if 6513 if (.not.btdb_get(rtdb,'nwpw:mc_group_size', 6514 > mt_int,1,mc_group_size)) then 6515 mc_group_size = 0 6516 end if 6517 6518 return 6519 end 6520 6521 6522 6523* *********************************** 6524* * * 6525* * control_mc_groups * 6526* * * 6527* *********************************** 6528 subroutine control_mc_groups(mc_group_start, 6529 > mc_group_end, 6530 > mc_group) 6531 implicit none 6532 integer mc_group_start(*) 6533 integer mc_group_end(*) 6534 integer mc_group(*) 6535 6536#include "bafdecls.fh" 6537#include "btdb.fh" 6538 6539* **** control_rtdb common block **** 6540 integer rtdb 6541 common / control_rtdb1 / rtdb 6542 6543 integer ng,ngs 6544 6545 if (.not.btdb_get(rtdb,'nwpw:mc_ngroups',mt_int,1,ng)) ng=1 6546 if (.not.btdb_get(rtdb,'nwpw:mc_group_size',mt_int,1,ngs)) ngs=1 6547 6548 if (.not.btdb_get(rtdb,'nwpw:mc_group_start', 6549 > mt_int,ng,mc_group_start)) then 6550 call icopy(ng,0,1,mc_group_start,1) 6551 end if 6552 if (.not.btdb_get(rtdb,'nwpw:mc_group_end', 6553 > mt_int,ng,mc_group_end)) then 6554 call icopy(ng,0,1,mc_group_end,1) 6555 end if 6556 if (.not.btdb_get(rtdb,'nwpw:mc_group', 6557 > mt_int,ngs,mc_group)) then 6558 call icopy(ngs,0,1,mc_group,1) 6559 end if 6560 6561 return 6562 end 6563 6564 6565* ************************************ 6566* * * 6567* * control_hess_model * 6568* * * 6569* ************************************ 6570 logical function control_hess_model() 6571 implicit none 6572 6573#include "control.fh" 6574 6575 control_hess_model = hess_model 6576 return 6577 end 6578 6579* ************************************ 6580* * * 6581* * control_hess_filename * 6582* * * 6583* ************************************ 6584 subroutine control_hess_filename(filehess) 6585 implicit none 6586 character*(*) filehess 6587 6588#include "control.fh" 6589#include "bafdecls.fh" 6590#include "btdb.fh" 6591#include "util.fh" 6592 6593* **** control_rtdb common block **** 6594 integer rtdb 6595 common / control_rtdb1 / rtdb 6596 6597 if (.not.btdb_cget(rtdb,'nwpw:hess_model:filename', 6598 > 1,filehess)) then 6599 call util_file_name('hess',.false.,.false.,filehess) 6600 end if 6601 return 6602 end 6603 6604 6605 6606* *********************************** 6607* * * 6608* * control_psp_semicore_small * 6609* * * 6610* *********************************** 6611 logical function control_psp_semicore_small() 6612 implicit none 6613 6614#include "bafdecls.fh" 6615#include "btdb.fh" 6616 6617* **** control_rtdb common block **** 6618 integer rtdb 6619 common / control_rtdb1 / rtdb 6620 6621 logical value 6622 6623 value = .true. 6624 if (.not.btdb_get(rtdb,'nwpw:psp:semicore_small',mt_log,1,value)) 6625 > value = .false. 6626 6627 control_psp_semicore_small = value 6628 return 6629 end 6630 6631 6632 6633 6634 6635 6636* *********************************** 6637* * * 6638* * control_fetch_confine * 6639* * * 6640* *********************************** 6641 subroutine control_fetch_confine(symbol,q1,rs,re) 6642 implicit none 6643 character*4 symbol 6644 real*8 q1,rs,re 6645 6646#include "bafdecls.fh" 6647#include "btdb.fh" 6648 6649* **** control_rtdb common block **** 6650 integer rtdb 6651 common / control_rtdb1 / rtdb 6652 6653 logical value 6654 integer l 6655 real*8 rrr(3) 6656 6657 l = index(symbol,' ') - 1 6658 if (.not.btdb_get(rtdb,'nwpw:confine:'//symbol(1:l), 6659 > mt_dbl,3,rrr)) then 6660 6661 rrr(1) = 1.0d0 6662 rrr(2) = 6.0d0 6663 rrr(3) = 14.0d0 6664 end if 6665 q1 = rrr(1) 6666 rs = rrr(2) 6667 re = rrr(3) 6668 return 6669 end 6670 6671 6672* *********************************** 6673* * * 6674* * control_runsocket * 6675* * * 6676* *********************************** 6677 logical function control_runsocket() 6678 implicit none 6679 6680#include "inp.fh" 6681 6682* **** local variables **** 6683 logical value 6684 character*40 stype 6685 6686 value = .false. 6687 call control_socket_type(stype) 6688 if (inp_compare(.false.,'ipi_client',stype)) value = .true. 6689 if (inp_compare(.false.,'unix',stype)) value = .true. 6690 6691 control_runsocket = value 6692 return 6693 end 6694 6695 6696* ************************************ 6697* * * 6698* * control_socket_type * 6699* * * 6700* ************************************ 6701 subroutine control_socket_type(stype) 6702 implicit none 6703 character*(*) stype 6704 6705#include "control.fh" 6706#include "btdb.fh" 6707 6708* **** control_rtdb common block **** 6709 integer rtdb 6710 common / control_rtdb1 / rtdb 6711 6712 if (.not.btdb_cget(rtdb,'nwpw:socket_type', 6713 > 1,stype)) then 6714 stype = "no socket" 6715 end if 6716 return 6717 end 6718 6719* ************************************ 6720* * * 6721* * control_socket_ip * 6722* * * 6723* ************************************ 6724 subroutine control_socket_ip(socket_ip) 6725 implicit none 6726 character*(*) socket_ip 6727 character*40 stype 6728 6729#include "control.fh" 6730#include "btdb.fh" 6731#include "inp.fh" 6732 6733* **** control_rtdb common block **** 6734 integer rtdb 6735 common / control_rtdb1 / rtdb 6736 6737 if (.not.btdb_cget(rtdb,'nwpw:socket_ip', 6738 > 1,socket_ip)) then 6739 call control_socket_type(stype) 6740 if (inp_compare(.false.,'unix',stype)) then 6741 socket_ip = "nwchem" 6742 else 6743 socket_ip = "127.0.0.1:31415" 6744 end if 6745 end if 6746 return 6747 end 6748 6749 6750 6751 6752 6753