1c 2c $Id$ 3c 4 5 6 subroutine cpmd_input(rtdb) 7 implicit none 8 integer rtdb 9#include "inp.fh" 10#include "bafdecls.fh" 11#include "rtdb.fh" 12#include "nwc_const.fh" 13#include "errquit.fh" 14c 15 logical value 16c 17 integer ind ! Index of matched directive 18 integer num_dirs ! No. of known directives 19 parameter (num_dirs = 61) 20 21 character*32 dirs(num_dirs) 22 character*255 test, id 23 24 data dirs / 'cell_name:', 25 > 'cell_name', 26 > 'geometry_optimize', 27 > 'formatted_filename:', 28 > 'formatted_filename ', 29 > 'input_wavefunction_filename:', 30 > 'input_wavefunction_filename ', 31 > 'output_wavefunction_filename:', 32 > 'output_wavefunction_filename', 33 > 'fake_mass:', 34 > 'fake_mass', 35 > 'time_step:', 36 > 'time_step', 37 > 'loop:', 38 > 'loop', 39 > 'scaling:', 40 > 'scaling', 41 > 'energy_cutoff:', 42 > 'energy_cutoff', 43 > 'wavefunction_cutoff:', 44 > 'wavefunction_cutoff', 45 > 'input_v_wavefunction_filename:', 46 > 'input_v_wavefunction_filename', 47 > 'output_v_wavefunction_filename:', 48 > 'output_v_wavefunction_filename', 49 > 'ewald_rcut:', 50 > 'ewald_rcut', 51 > 'ewald_ncut:', 52 > 'ewald_ncut', 53 > 'xyz_filename:', 54 > 'xyz_filename', 55 > 'exchange_correlation:', 56 > 'exchange_correlation', 57 > 'xc', 58 > 'fractional_coordinates', 59 > 'Nose-Hoover:', 60 > 'Nose-Hoover', 61 > 'emotion_filename:', 62 > 'emotion_filename', 63 > 'eigmotion_filename:', 64 > 'eigmotion_filename', 65 > 'hmotion_filename:', 66 > 'hmotion_filename', 67 > 'omotion_filename:', 68 > 'omotion_filename', 69 > 'ion_motion_filename:', 70 > 'ion_motion_filename', 71 > 'mulliken', 72 > 'energy', 73 > 'SA_decay', 74 > 'Fei', 75 > 'Fei_quench', 76 > 'pressure', 77 > 'cutoff', 78 > 'initial_velocities', 79 > 'reset_Temperature', 80 > 'rotation', 81 > 'translation', 82 > 'dipole_motion', 83 > 'temperature', 84 > 'end'/ 85 86 character*50 cell_name 87 character*50 input_wavefunction_filename 88 character*50 output_wavefunction_filename 89 character*50 input_v_wavefunction_filename 90 character*50 output_v_wavefunction_filename 91 character*50 xyz_filename 92 character*50 emotion_filename 93 character*50 eigmotion_filename 94 character*50 hmotion_filename 95 character*50 omotion_filename 96 character*50 ion_motion_filename 97 character*50 exchange_correlation,zone_name 98 logical geometry_optimize,frac_coord,mulliken,mulliken_kawai 99 double precision fake_mass,time_step,rcut 100 integer loop(2),npsp,ncut,mchain,nchain 101 double precision scaling(2),ecut,wcut 102 logical nose,nosers,move 103 double precision Pe,Te,fe 104 double precision Pr,Tr 105 integer index_start1(2),n1,j,jstart,jlast,jstride 106 107 logical nwpw_parse_boolean 108 external nwpw_parse_boolean 109 110 111* ***** initializations **** 112c call cpmd_input_default(rtdb) 113 npsp = 0 114 115 116 10 if (.not. inp_read()) 117 > call errquit( 118 > 'cpmd_input: inp_read failed', 0, INPUT_ERR) 119 if (.not. inp_a(test)) 120 > call errquit( 121 > 'cpmd_input: failed to read keyword', 0, INPUT_ERR) 122 if (.not. inp_match(num_dirs, .false., test, dirs, ind)) 123 > call errquit( 124 > 'cpmd_input: unknown directive', 0, INPUT_ERR) 125 126 127 goto ( 100,100, 200, 300,300, 400,400, 500,500, 600,600, 128 > 700,700, 800,800, 900,900, 1000,1000, 1100,1100, 129 > 1200,1200, 1300,1300, 1400,1400, 1500,1500, 1600,1600, 130 > 1700,1700,1700, 1800, 1900,1900, 2000,2000, 2100,2100, 131 > 2200,2200, 2300,2300, 2400,2400, 132 > 2500,2600,2700,2800,2900,3000,1150,3100,3200, 133 > 3300,3400,3500,3600, 134 > 9999) ind 135 call errquit( 136 > 'psp_formatter_input: unimplemented directive', ind, 137 & INPUT_ERR) 138 139c 140c cell_name 141c 142 100 if (.not. inp_a(cell_name)) 143 > call errquit( 144 > 'cpmd_input: failed to read cell_name', 0, INPUT_ERR) 145 ind = index(cell_name,' ') - 1 146 value = rtdb_cput(rtdb,'cpmd:cell_name',1,cell_name(1:ind)) 147 if (.not.value) 148 > call errquit('cpmd_input: writing cell_name', 100, RTDB_ERR) 149 goto 10 150c 151c geometry_optimize 152c 153 200 geometry_optimize = .true. 154 value = rtdb_put(rtdb,'cpmd:geometry_optimize',mt_log,1, 155 > geometry_optimize) 156 if (.not.value) 157 > call errquit('cpmd_input: writing ', 200, RTDB_ERR) 158 goto 10 159c 160c formatted_psp_filename 161c 162 300 if (.not. inp_a(test)) 163 > call errquit( 164 > 'cpmd_input: failed to read psp_filename', 0, INPUT_ERR) 165 npsp = npsp + 1 166 id = 'cpmd:psp'//CHAR(npsp) 167 ind = index(test,' ') - 1 168 if (.not. rtdb_cput(rtdb,id, 169 > 1,test(1:ind))) 170 > call errquit( 171 > 'cpmd_input: rtdb_cput failed', 0, RTDB_ERR) 172 value = rtdb_put(rtdb,'cpmd:npsp', mt_int,1,npsp) 173 if (.not.value) 174 > call errquit('cpmd_input: writing ', 300, RTDB_ERR) 175 goto 10 176c 177c input_wavefunction_filename 178c 179 400 if (.not. inp_a(input_wavefunction_filename)) 180 > call errquit( 181 > 'cpmd_input: failed to read input_wavefunctions', 0, 182 & INPUT_ERR) 183 ind = index(input_wavefunction_filename,' ') - 1 184 value = rtdb_cput(rtdb,'cpmd:input_wavefunction_filename', 185 > 1,input_wavefunction_filename(1:ind)) 186 if (.not.value) 187 > call errquit('cpmd_input: writing ', 400, RTDB_ERR) 188 goto 10 189c 190c output_wavefunction_filename 191c 192 500 if (.not. inp_a(output_wavefunction_filename)) 193 > call errquit( 194 > 'cpmd_input: failed to read output_wavefunction', 0, 195 & INPUT_ERR) 196 ind = index(output_wavefunction_filename,' ') - 1 197 value = rtdb_cput(rtdb,'cpmd:output_wavefunction_filename', 198 > 1,output_wavefunction_filename(1:ind)) 199 if (.not.value) 200 > call errquit('cpmd_input: writing ', 500, RTDB_ERR) 201 goto 10 202c 203c fake_mass 204c 205 600 if (.not. inp_f(fake_mass)) 206 > call errquit( 207 > 'cpmd_input: failed to read fake_mass', 0, INPUT_ERR) 208 value = rtdb_put(rtdb,'cpmd:fake_mass',mt_dbl,1,fake_mass) 209 if (.not.value) 210 > call errquit('cpmd_input: writing ', 600, RTDB_ERR) 211 goto 10 212c 213c time_step 214c 215 700 if (.not. inp_f(time_step)) 216 > call errquit( 217 > 'cpmd_input: failed to read time_step', 0, INPUT_ERR) 218 value = rtdb_put(rtdb,'cpmd:time_step',mt_dbl,1,time_step) 219 if (.not.value) 220 > call errquit('cpmd_input: writing ', 700, RTDB_ERR) 221 goto 10 222c 223c loop 224c 225 800 if (.not. inp_i(loop(1))) 226 > call errquit( 227 > 'cpmd_input: failed to read loop', 0, INPUT_ERR) 228 if (.not. inp_i(loop(2))) 229 > call errquit( 230 > 'cpmd_input: failed to read loop', 0, INPUT_ERR) 231 value = rtdb_put(rtdb,'cpmd:loop',mt_int,2,loop) 232 if (.not.value) 233 > call errquit('cpmd_input: writing ', 800, RTDB_ERR) 234 goto 10 235c 236c scaling 237c 238 900 if (.not. inp_f(scaling(1))) 239 > call errquit( 240 > 'cpmd_input: failed to read scaling', 0, INPUT_ERR) 241 if (.not. inp_f(scaling(2))) 242 > call errquit( 243 > 'cpmd_input: failed to read scaling', 0, INPUT_ERR) 244 value = rtdb_put(rtdb,'cpmd:scaling',mt_dbl,2,scaling) 245 246 !*** try reading in atom indexes **** 247 if (.not.BA_push_get(mt_int,nw_max_atom,'indx_start1', 248 > index_start1(2),index_start1(1))) 249 > call errquit( 250 > 'cpmd_input:failed allocating index_start1',0,MA_ERR) 251 n1 = 0 252 do while (inp_irange(jstart,jlast,jstride)) 253 do j=jstart,jlast,jstride 254 int_mb(index_start1(1)+n1) = j 255 n1 = n1+1 256 end do 257 end do 258 value = value.and.rtdb_put(rtdb,'nwpw:scaling_natms', 259 > mt_int,1,n1) 260 if (n1.gt.0) then 261 value = value.and.rtdb_put(rtdb,'nwpw:scaling_atoms', 262 > mt_int,n1,int_mb(index_start1(1))) 263 end if 264 265 if (.not.BA_pop_stack(index_start1(2))) call errquit( 266 > 'cpmd_input:failed deallocating index_start1',0,MA_ERR) 267 268 if (.not.value) 269 > call errquit('cpmd_input: writing ', 900, RTDB_ERR) 270 goto 10 271c 272c energy_cutoff 273c 274 1000 if (.not. inp_f(ecut)) 275 > call errquit( 276 > 'cpmd_input: failed to read ecut', 0, INPUT_ERR) 277 value = rtdb_put(rtdb,'cpmd:ecut',mt_dbl,1,ecut) 278 if (.not.value) 279 > call errquit('cpmd_input: writing ', 1000, RTDB_ERR) 280 goto 10 281c 282c wavefunction_cutoff 283c 284 1100 if (.not. inp_f(wcut)) 285 > call errquit( 286 > 'cpmd_input: failed to read wcut', 0, INPUT_ERR) 287 value = rtdb_put(rtdb,'cpmd:wcut',mt_dbl,1,wcut) 288 if (.not.value) 289 > call errquit('cpmd_input: writing ', 1100, RTDB_ERR) 290 goto 10 291 292c 293c cutoff 294c 295 1150 if (.not. inp_f(wcut)) 296 > call errquit( 297 > 'cpmd_input: failed to read wcut', 0, INPUT_ERR) 298 value = rtdb_put(rtdb,'cpmd:wcut',mt_dbl,1,wcut) 299 value = rtdb_put(rtdb,'cpmd:ecut',mt_dbl,1,2.0d0*wcut) 300 if (.not.value) 301 > call errquit('cpmd_input: writing ', 1150, RTDB_ERR) 302 goto 10 303 304c 305c input_v_wavefunction_filename 306c 307 1200 if (.not. inp_a(input_v_wavefunction_filename)) 308 > call errquit( 309 > 'cpmd_input: failed to read input_v_wavefunction', 0, 310 & INPUT_ERR) 311 ind = index(input_v_wavefunction_filename,' ') - 1 312 value = rtdb_cput(rtdb,'cpmd:input_v_wavefunction_filename', 313 > 1,input_v_wavefunction_filename(1:ind)) 314 if (.not.value) 315 > call errquit('cpmd_input: writing ', 1200, RTDB_ERR) 316 317 goto 10 318c 319c output_v_wavefunction_filename 320c 321 1300 if (.not. inp_a(output_v_wavefunction_filename)) 322 > call errquit( 323 > 'cpmd_input: failed to read output_v_wavefunction', 0, 324 & INPUT_ERR) 325 ind = index(output_v_wavefunction_filename,' ') - 1 326 value = rtdb_cput(rtdb,'cpmd:output_v_wavefunction_filename', 327 > 1,output_v_wavefunction_filename(1:ind)) 328 if (.not.value) 329 > call errquit('cpmd_input: writing ', 1300, RTDB_ERR) 330 goto 10 331c 332c ewald_rcut 333c 334 1400 if (.not. inp_f(rcut)) 335 > call errquit( 336 > 'cpmd_input: failed to read rcut', 0, INPUT_ERR) 337 value = rtdb_put(rtdb,'cpmd:rcut',mt_dbl,1,rcut) 338 if (.not.value) 339 > call errquit('cpmd_input: writing ', 1400, RTDB_ERR) 340 goto 10 341c 342c ewald_ncut 343c 344 1500 if (.not. inp_i(ncut)) 345 > call errquit( 346 > 'cpmd_input: failed to read ncut', 0, INPUT_ERR) 347 value = rtdb_put(rtdb,'cpmd:ncut',mt_int,1,ncut) 348 if (.not.value) 349 > call errquit('cpmd_input: writing ', 1500, RTDB_ERR) 350 goto 10 351c 352c xyz_filename 353c 354 1600 if (.not. inp_a(xyz_filename)) 355 > call errquit( 356 > 'cpmd_input: failed to read xyz_filename', 0, INPUT_ERR) 357 ind = index(xyz_filename,' ') - 1 358 value = rtdb_cput(rtdb,'cpmd:xyz_filename', 359 > 1,xyz_filename(1:ind)) 360 if (.not.value) 361 > call errquit('cpmd_input: writing ', 1600, RTDB_ERR) 362 goto 10 363c 364c exchange_correlation 365c 366 1700 if (.not. inp_a(exchange_correlation)) 367 > call errquit( 368 > 'cpmd_input: failed to read exchange_correlation', 0, 369 & INPUT_ERR) 370 ind = index(exchange_correlation,' ') - 1 371 value = rtdb_cput(rtdb,'cpmd:exchange_correlation', 372 > 1,exchange_correlation(1:ind)) 373 if (.not.value) 374 > call errquit('cpmd_input: writing ', 1700, RTDB_ERR) 375 goto 10 376c 377c fractional_coordinates 378c 379 1800 frac_coord = .true. 380 value = rtdb_put(rtdb,'cpmd:fractional_coordinates',mt_log,1, 381 > frac_coord) 382 if (.not.value) 383 > call errquit('cpmd_input: writing ', 1800, RTDB_ERR) 384 goto 10 385c 386c Nose-Hoover 387c 388 1900 nose = .true. 389 nosers = .true. 390 Pe = 1200.00d0 391 Pr = 1200.00d0 392 Te = 298.15d0 393 Tr = 298.15d0 394 mchain = 1 395 nchain = 1 396 if (.not. inp_f(fe)) goto 1901 397 Pe = fe 398 if (.not. inp_f(fe)) goto 1901 399 Te = fe 400 if (.not. inp_f(fe)) goto 1901 401 Pr = fe 402 if (.not. inp_f(fe)) goto 1901 403 Tr = fe 404 if (.not.inp_a(zone_name)) goto 1901 405 if (inp_compare(.false.,zone_name,'restart')) nosers=.true. 406 if (inp_compare(.false.,zone_name,'start')) nosers=.false. 407 if (inp_compare(.false.,zone_name,'clear')) nosers=.false. 408 if (inp_compare(.false.,zone_name,'initialize')) nosers=.false. 409 if (.not.inp_i(mchain)) goto 1901 410 if (.not.inp_i(nchain)) goto 1901 411 1901 value = rtdb_put(rtdb,'cpmd:nose',mt_log,1,nose) 412 value = value.and. 413 > rtdb_put(rtdb,'cpmd:nose_restart',mt_log,1,nosers) 414 value = value.and.rtdb_put(rtdb,'cpmd:Pe',mt_dbl,1,Pe) 415 value = value.and.rtdb_put(rtdb,'cpmd:Te',mt_dbl,1,Te) 416 value = value.and.rtdb_put(rtdb,'cpmd:Pr',mt_dbl,1,Pr) 417 value = value.and.rtdb_put(rtdb,'cpmd:Tr',mt_dbl,1,Tr) 418 value = value.and.rtdb_put(rtdb,'cpmd:Mchain',mt_int,1,mchain) 419 value = value.and.rtdb_put(rtdb,'cpmd:Nchain',mt_int,1,nchain) 420 if (.not.value) 421 > call errquit('cpmd_input: writing ', 1900, RTDB_ERR) 422 goto 10 423c 424c emotion_filename 425c 426 2000 if (.not. inp_a(emotion_filename)) 427 > call errquit( 428 > 'cpmd_input: failed to read emotion_filename', 0, 429 & INPUT_ERR) 430 ind = index(emotion_filename,' ') - 1 431 value = rtdb_cput(rtdb,'cpmd:emotion_filename', 432 > 1,emotion_filename(1:ind)) 433 if (.not.value) 434 > call errquit('cpmd_input: writing ', 2000, RTDB_ERR) 435 goto 10 436c 437c eigmotion_filename 438c 439 2100 if (.not. inp_a(eigmotion_filename)) 440 > call errquit( 441 > 'cpmd_input: failed to read eigmotion_filename', 0, 442 & INPUT_ERR) 443 ind = index(eigmotion_filename,' ') - 1 444 value = rtdb_cput(rtdb,'cpmd:eigmotion_filename', 445 > 1,eigmotion_filename(1:ind)) 446 if (.not.value) 447 > call errquit('cpmd_input: writing ', 2100, RTDB_ERR) 448 goto 10 449c 450c hmotion_filename 451c 452 2200 if (.not. inp_a(hmotion_filename)) 453 > call errquit( 454 > 'cpmd_input: failed to read hmotion_filename', 0, 455 & INPUT_ERR) 456 ind = index(hmotion_filename,' ') - 1 457 value = rtdb_cput(rtdb,'cpmd:hmotion_filename', 458 > 1,hmotion_filename(1:ind)) 459 if (.not.value) 460 > call errquit('cpmd_input: writing ', 2200, RTDB_ERR) 461 goto 10 462c 463c omotion_filename 464c 465 2300 if (.not. inp_a(omotion_filename)) 466 > call errquit( 467 > 'cpmd_input: failed to read omotion_filename', 0, 468 & INPUT_ERR) 469 ind = index(omotion_filename,' ') - 1 470 value = rtdb_cput(rtdb,'cpmd:omotion_filename', 471 > 1,omotion_filename(1:ind)) 472 if (.not.value) 473 > call errquit('cpmd_input: writing ', 2300, RTDB_ERR) 474 goto 10 475c 476c ion_motion_filename 477c 478 2400 if (.not. inp_a(ion_motion_filename)) 479 > call errquit( 480 > 'cpmd_input: failed to read ion_motion_filename', 0, 481 & INPUT_ERR) 482 ind = index(ion_motion_filename,' ') - 1 483 value = rtdb_cput(rtdb,'cpmd:ion_motion_filename', 484 > 1,ion_motion_filename(1:ind)) 485 if (.not.value) 486 > call errquit('cpmd_input: writing ', 2400, RTDB_ERR) 487 goto 10 488c 489c Mulliken 490c 491 2500 if (.not.inp_a(zone_name)) zone_name = 'lcao' 492 mulliken = .true. 493 mulliken_kawai = .false. 494 ind = index(zone_name,' ') - 1 495 if (inp_compare(.false.,zone_name,'kawai')) mulliken_kawai=.true. 496 value = rtdb_put(rtdb,'cpmd:mulliken',mt_log,1, 497 > mulliken) 498 > .and.rtdb_put(rtdb,'nwpw:mulliken_kawai',mt_log,1, 499 > mulliken_kawai) 500 if (.not.value) 501 > call errquit( 502 > 'cpmd_input: error writing to rtdb', 2500, RTDB_ERR) 503 goto 10 504 505c 506c energy 507c 508 2600 nose = .false. 509 value = rtdb_put(rtdb,'cpmd:nose',mt_log,1,nose) 510 if (.not.value) 511 > call errquit('cpmd_input: writing ', 2600, RTDB_ERR) 512 goto 10 513c 514c SA_decay 515c 516 2700 if (inp_f(scaling(1))) then 517 scaling(2) = scaling(1) 518 if (.not.inp_f(scaling(2))) scaling(2) = scaling(1) 519 value = rtdb_put(rtdb,'cpmd:sa_decay',mt_dbl,2,scaling) 520 else 521 value = rtdb_delete(rtdb,'cpmd:sa_decay') 522 end if 523 if (.not.value) 524 > call errquit('cpmd_input: writing sa_decay', 2700, RTDB_ERR) 525 goto 10 526c 527c Fei 528c 529 2800 if (inp_a(ion_motion_filename)) then 530 ind = index(ion_motion_filename,' ') - 1 531 value = rtdb_cput(rtdb,'cpmd:fei_filename', 532 > 1,ion_motion_filename(1:ind)) 533 value = rtdb_cput(rtdb,'nwpw:fei_filename', 534 > 1,ion_motion_filename(1:ind)) 535 if (.not.value) 536 > call errquit('cpmd_input: writing ', 2800, RTDB_ERR) 537 end if 538 value = .true. 539 if (.not.rtdb_put(rtdb,'cpmd:fei',mt_log,1,value)) 540 > call errquit('cpmd_input: writing fei ', 2800, RTDB_ERR) 541 if (.not.rtdb_put(rtdb,'nwpw:fei',mt_log,1,value)) 542 > call errquit('cpmd_input: writing fei ', 2800, RTDB_ERR) 543 goto 10 544c 545c Fei Quench 546c 547 2900 value = .true. 548 if (.not.rtdb_put(rtdb,'cpmd:fei_quench',mt_log,1,value)) 549 > call errquit('cpmd_input: writing fei quench ', 2900, RTDB_ERR) 550 goto 10 551 552c 553c pressure 554c 555 3000 if (.not.inp_a(zone_name)) zone_name = 'on' 556 value = .true. 557 value = nwpw_parse_boolean(zone_name,value) 558c if (inp_compare(.false.,zone_name,'on')) value=.true. 559c if (inp_compare(.false.,zone_name,'.true.')) value=.true. 560c if (inp_compare(.false.,zone_name,'off')) value=.false. 561c if (inp_compare(.false.,zone_name,'.false.')) value=.false. 562 if (.not.rtdb_put(rtdb,'cpmd:pressure',mt_log,1,value)) 563 > call errquit('cpmd_input: writing pressure ', 3000, RTDB_ERR) 564 goto 10 565 566c 567c initial_velocities temperature seed 568c rtdb = cpmd:init_velocities_temperature 569c 570 3100 if (inp_f(scaling(1))) then !*** in units of K *** 571 if (.not.inp_i(n1)) n1 = 494 572 573 if (.not.rtdb_put(rtdb,'nwpw:init_velocities_temperature', 574 > mt_dbl,1,scaling)) 575 > call errquit('cpmd_input: initial temperature',3100,RTDB_ERR) 576 577 if (.not.rtdb_put(rtdb,'nwpw:init_velocities_seed', 578 > mt_int,1,n1)) 579 > call errquit('cpmd_input: initial temperature',3100,RTDB_ERR) 580 581 if (.not.rtdb_put(rtdb,'nwpw:init_velocities', 582 > mt_log,1,.true.)) 583 > call errquit('cpmd_input: initial velocities',3100,RTDB_ERR) 584 585 else 586 value = rtdb_delete(rtdb,'nwpw:init_velocities') 587 end if 588 goto 10 589 590 3200 value = rtdb_delete(rtdb,"nwpw:ke_total") 591 value = rtdb_delete(rtdb,"nwpw:kg_total") 592 value = rtdb_delete(rtdb,"nwpw:ke_count") 593 goto 10 594 595 596c 597c rotation 598c 599 3300 if (.not.inp_a(zone_name)) goto 3301 600 move = nwpw_parse_boolean(zone_name,.false.) 601 if (.not.rtdb_put(rtdb,'nwpw:rotation',mt_log,1,move)) 602 > call errquit('cpmd_input: error writing to rtdb',3300,RTDB_ERR) 603 604 3301 continue 605 goto 10 606 607c 608c translation 609c 610 3400 if (.not.inp_a(zone_name)) goto 3401 611 move = nwpw_parse_boolean(zone_name,.true.) 612 if (.not.rtdb_put(rtdb,'cgsd:allow_translation',mt_log,1,move)) 613 > call errquit('cpmd_input: error writing to rtdb',3400,RTDB_ERR) 614 if (.not.rtdb_put(rtdb,'band:allow_translation',mt_log,1,move)) 615 > call errquit('cpmd_input: error writing to rtdb',3400,RTDB_ERR) 616 617 3401 continue 618 goto 10 619 620c 621c dipole_motion 622c 623 3500 if (inp_a(ion_motion_filename)) then 624 ind = index(ion_motion_filename,' ') - 1 625 if (.not.rtdb_cput(rtdb,'nwpw:dipole_motion_filename', 626 > 1,ion_motion_filename(1:ind))) 627 > call errquit('nwpw_input: writing filename',3500,RTDB_ERR) 628 end if 629 if (.not.rtdb_put(rtdb,'nwpw:dipole_motion',mt_log,1,value)) 630 > call errquit('cpmd_input: writing dipole_motion',3500,RTDB_ERR) 631 goto 10 632 633 634c 635c temperature Tion Pion Telc Pion Pelc start/restart nchain mchain 636c - a reordering of the parameters for nose-hoover 637c 638 3600 nose = .true. 639 nosers = .true. 640 Tr = 298.15d0 641 Pr = 1200.00d0 642 Te = 298.15d0 643 Pe = 1200.00d0 644 mchain = 1 645 nchain = 1 646 if (.not. inp_f(fe)) goto 3601 647 Tr = fe 648 Te = fe 649 if (.not. inp_f(fe)) goto 3601 650 Pr = fe 651 Pe = fe 652 if (.not. inp_f(fe)) goto 3601 653 Te = fe 654 if (.not. inp_f(fe)) goto 3601 655 Pe = fe 656 3601 if (.not.inp_a(zone_name)) goto 3602 657 if (inp_compare(.false.,zone_name,'restart')) nosers=.true. 658 if (inp_compare(.false.,zone_name,'start')) nosers=.false. 659 if (inp_compare(.false.,zone_name,'clear')) nosers=.false. 660 if (inp_compare(.false.,zone_name,'initialize')) nosers=.false. 661 3602 if (.not.inp_i(nchain)) goto 3603 662 mchain = nchain 663 if (.not.inp_i(mchain)) goto 3603 664 665 3603 value = rtdb_put(rtdb,'cpmd:nose',mt_log,1,nose) 666 value = value.and. 667 > rtdb_put(rtdb,'cpmd:nose_restart',mt_log,1,nosers) 668 value = value.and.rtdb_put(rtdb,'cpmd:Pe',mt_dbl,1,Pe) 669 value = value.and.rtdb_put(rtdb,'cpmd:Te',mt_dbl,1,Te) 670 value = value.and.rtdb_put(rtdb,'cpmd:Pr',mt_dbl,1,Pr) 671 value = value.and.rtdb_put(rtdb,'cpmd:Tr',mt_dbl,1,Tr) 672 value = value.and.rtdb_put(rtdb,'cpmd:Mchain',mt_int,1,mchain) 673 value = value.and.rtdb_put(rtdb,'cpmd:Nchain',mt_int,1,nchain) 674 if (.not.value) 675 > call errquit('cpmd_input: writing ', 3600, RTDB_ERR) 676 goto 10 677 678 679 680* ***** add wavefunction to rtdb **** 681 9999 continue 682 return 683 end 684 685 subroutine cpmd_input_default(rtdb) 686 implicit none 687#include "errquit.fh" 688 integer rtdb 689#include "inp.fh" 690#include "bafdecls.fh" 691#include "rtdb.fh" 692 693* ***** local variables ***** 694 logical value 695 integer ind 696 697 character*50 cell_name 698 character*50 input_wavefunction_filename 699 character*50 output_wavefunction_filename 700 character*50 input_v_wavefunction_filename 701 character*50 output_v_wavefunction_filename 702 character*50 xyz_filename 703 character*50 emotion_filename 704 character*50 eigmotion_filename 705 character*50 hmotion_filename 706 character*50 omotion_filename 707 character*50 ion_motion_filename 708 character*50 exchange_correlation 709 logical geometry_optimize,frac_coord,mulliken 710 double precision fake_mass,time_step,rcut 711 integer loop(2),npsp,ncut 712 double precision scaling(2),ecut,wcut 713 logical nose 714 double precision Pe,Te 715 double precision Pr,Tr 716 717 718* **** don't set defaults if they already exist **** 719 value = rtdb_get(rtdb,'cpmd:ncut',mt_int,1,ncut) 720 if (value) return 721 722 723* ***** initializations **** 724 cell_name = 'cell_default' 725 726 call util_file_prefix('movecs',input_wavefunction_filename) 727 call util_file_prefix('movecs',output_wavefunction_filename) 728c input_wavefunction_filename = ' ' 729c output_wavefunction_filename = ' ' 730 731 call util_file_prefix('vmovecs',input_v_wavefunction_filename) 732 call util_file_prefix('vmovecs',output_v_wavefunction_filename) 733c input_v_wavefunction_filename = ' ' 734c output_v_wavefunction_filename = ' ' 735 736 !exchange_correlation = 'vosko' 737 geometry_optimize = .true. 738 frac_coord = .false. 739 mulliken = .false. 740 fake_mass = 1000.0d0 741 time_step = 5.0d0 742 loop(1) = 10 743 loop(2) = 1 744 scaling(1) = 1.0d0 745 scaling(2) = 1.0d0 746 ecut=9000.0d0 747 wcut=9000.0d0 748 rcut=0.0d0 749 ncut=0 750 npsp = 0 751 752 call util_file_prefix('xyz',xyz_filename) 753 call util_file_prefix('emotion',emotion_filename) 754 call util_file_prefix('eigmotion',eigmotion_filename) 755 call util_file_prefix('hmotion',hmotion_filename) 756 call util_file_prefix('omotion',omotion_filename) 757 call util_file_prefix('ion_motion',ion_motion_filename) 758c xyz_filename = 'XYZ' 759c emotion_filename = 'EMOTION' 760c eigmotion_filename = 'EIGMOTION' 761c hmotion_filename = 'HMOTION' 762c omotion_filename = 'OMOTION' 763c ion_motion_filename = 'MOTION' 764 765 nose = .false. 766 Pe = 1200.00d0 767 Pr = 1200.00d0 768 Te = 298.15d0 769 Tr = 298.15d0 770 771 ind = index(cell_name,' ') - 1 772 value = rtdb_cput(rtdb,'cpmd:cell_name',1,cell_name(1:ind)) 773 774 ind = index(input_wavefunction_filename,' ') - 1 775 value = value.and. 776 > rtdb_cput(rtdb,'cpmd:input_wavefunction_filename', 777 > 1,input_wavefunction_filename(1:ind)) 778 779 ind = index(output_wavefunction_filename,' ') - 1 780 value = value.and. 781 > rtdb_cput(rtdb,'cpmd:output_wavefunction_filename', 782 > 1,output_wavefunction_filename(1:ind)) 783 784 ind = index(input_v_wavefunction_filename,' ') - 1 785 value = value.and. 786 > rtdb_cput(rtdb,'cpmd:input_v_wavefunction_filename', 787 > 1,input_v_wavefunction_filename(1:ind)) 788 789 ind = index(output_v_wavefunction_filename,' ') - 1 790 value = value.and. 791 > rtdb_cput(rtdb,'cpmd:output_v_wavefunction_filename', 792 > 1,output_v_wavefunction_filename(1:ind)) 793 794 ind = index(xyz_filename,' ') - 1 795 value = value.and. 796 > rtdb_cput(rtdb,'cpmd:xyz_filename', 797 > 1,xyz_filename(1:ind)) 798 799c ind = index(exchange_correlation,' ') - 1 800c value = value.and. 801c > rtdb_cput(rtdb,'cpmd:exchange_correlation', 802c > 1,exchange_correlation(1:ind)) 803 804 ind = index(emotion_filename,' ') - 1 805 value = value.and. 806 > rtdb_cput(rtdb,'cpmd:emotion_filename', 807 > 1,emotion_filename(1:ind)) 808 ind = index(eigmotion_filename,' ') - 1 809 value = value.and. 810 > rtdb_cput(rtdb,'cpmd:eigmotion_filename', 811 > 1,eigmotion_filename(1:ind)) 812 ind = index(hmotion_filename,' ') - 1 813 value = value.and. 814 > rtdb_cput(rtdb,'cpmd:hmotion_filename', 815 > 1,hmotion_filename(1:ind)) 816 ind = index(omotion_filename,' ') - 1 817 value = value.and. 818 > rtdb_cput(rtdb,'cpmd:omotion_filename', 819 > 1,omotion_filename(1:ind)) 820 ind = index(ion_motion_filename,' ') - 1 821 value = value.and. 822 > rtdb_cput(rtdb,'cpmd:ion_motion_filename', 823 > 1,ion_motion_filename(1:ind)) 824 825 value = value.and. 826 > rtdb_put(rtdb,'cpmd:geometry_optimize',mt_log,1, 827 > geometry_optimize) 828 value = value.and. 829 > rtdb_put(rtdb,'cpmd:fractional_coordinates',mt_log,1, 830 > frac_coord) 831 value = value.and. 832 > rtdb_put(rtdb,'cpmd:mulliken',mt_log,1, 833 > mulliken) 834 835 value = value.and. 836 > rtdb_put(rtdb,'cpmd:npsp', mt_int,1,npsp) 837 value = value.and. 838 > rtdb_put(rtdb,'cpmd:fake_mass',mt_dbl,1,fake_mass) 839 value = value.and. 840 > rtdb_put(rtdb,'cpmd:time_step',mt_dbl,1,time_step) 841 value = value.and. 842 > rtdb_put(rtdb,'cpmd:loop',mt_int,2,loop) 843 value = value.and. 844 > rtdb_put(rtdb,'cpmd:scaling',mt_dbl,2,scaling) 845 value = value.and. 846 > rtdb_put(rtdb,'cpmd:ecut',mt_dbl,1,ecut) 847 value = value.and. 848 > rtdb_put(rtdb,'cpmd:wcut',mt_dbl,1,wcut) 849 value = value.and. 850 > rtdb_put(rtdb,'cpmd:rcut',mt_dbl,1,rcut) 851 value = value.and. 852 > rtdb_put(rtdb,'cpmd:ncut',mt_int,1,ncut) 853 854 value = value.and. 855 > rtdb_put(rtdb,'cpmd:nose',mt_log,1,nose) 856 value = value.and. 857 > rtdb_put(rtdb,'cpmd:Pe',mt_dbl,1,Pe) 858 value = value.and. 859 > rtdb_put(rtdb,'cpmd:Te',mt_dbl,1,Te) 860 value = value.and. 861 > rtdb_put(rtdb,'cpmd:Pr',mt_dbl,1,Pr) 862 value = value.and. 863 > rtdb_put(rtdb,'cpmd:Tr',mt_dbl,1,Tr) 864 865* ***** Error writing to RTDB ***** 866 if (.not.value) 867 > call errquit( 868 > 'cpmd_default: error writing to rtdb', 0, RTDB_ERR) 869 870 return 871 end 872 873