1!--------------------------------------------------------------------------------------------------! 2! CP2K: A general program to perform molecular dynamics simulations ! 3! Copyright (C) 2000 - 2020 CP2K developers group ! 4!--------------------------------------------------------------------------------------------------! 5 6! ************************************************************************************************** 7!> \brief I/O subroutines for pint_env 8!> \author Lukasz Walewski 9!> \date 2009-06-04 10! ************************************************************************************************** 11MODULE pint_io 12 13 USE cell_types, ONLY: cell_type 14 USE cp_log_handling, ONLY: cp_get_default_logger,& 15 cp_logger_get_default_unit_nr,& 16 cp_logger_type 17 USE cp_output_handling, ONLY: cp_p_file,& 18 cp_print_key_finished_output,& 19 cp_print_key_should_output,& 20 cp_print_key_unit_nr 21 USE cp_subsys_types, ONLY: cp_subsys_get,& 22 cp_subsys_type 23 USE cp_units, ONLY: cp_unit_from_cp2k 24 USE f77_interface, ONLY: f_env_add_defaults,& 25 f_env_rm_defaults,& 26 f_env_type 27 USE force_env_types, ONLY: force_env_get 28 USE input_constants, ONLY: dump_atomic,& 29 dump_dcd,& 30 dump_dcd_aligned_cell,& 31 dump_xmol 32 USE input_section_types, ONLY: section_vals_get_subs_vals,& 33 section_vals_type,& 34 section_vals_val_get 35 USE kinds, ONLY: default_string_length,& 36 dp 37 USE machine, ONLY: m_flush 38 USE particle_list_types, ONLY: particle_list_type 39 USE particle_methods, ONLY: write_particle_coordinates 40 USE pint_public, ONLY: pint_com_pos 41 USE pint_transformations, ONLY: pint_u2x 42 USE pint_types, ONLY: e_conserved_id,& 43 e_kin_thermo_id,& 44 e_kin_virial_id,& 45 e_potential_id,& 46 pint_env_type 47#include "../base/base_uses.f90" 48 49 IMPLICIT NONE 50 51 PRIVATE 52 53 LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .TRUE. 54 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'pint_io' 55 56 PUBLIC :: pint_write_line 57 PUBLIC :: pint_write_centroids 58 PUBLIC :: pint_write_trajectory 59 PUBLIC :: pint_write_com 60 PUBLIC :: pint_write_ener 61 PUBLIC :: pint_write_action 62 PUBLIC :: pint_write_step_info 63 PUBLIC :: pint_write_rgyr 64 65CONTAINS 66 67! *************************************************************************** 68!> \brief Writes out a line of text to the default output unit. 69!> \param line ... 70!> \date 2009-07-10 71!> \author Lukasz Walewski 72! ************************************************************************************************** 73 SUBROUTINE pint_write_line(line) 74 75 CHARACTER(len=*), INTENT(IN) :: line 76 77 CHARACTER(len=*), PARAMETER :: routineN = 'pint_write_line', & 78 routineP = moduleN//':'//routineN 79 80 CHARACTER(len=default_string_length) :: my_label 81 INTEGER :: unit_nr 82 TYPE(cp_logger_type), POINTER :: logger 83 84 NULLIFY (logger) 85 logger => cp_get_default_logger() 86 my_label = "PINT|" 87 88 IF (logger%para_env%ionode) THEN 89 unit_nr = cp_logger_get_default_unit_nr(logger) 90 WRITE (unit_nr, '(T2,A)') TRIM(my_label)//" "//TRIM(line) 91 END IF 92 93 RETURN 94 END SUBROUTINE pint_write_line 95 96! *************************************************************************** 97!> \brief Write out the trajectory of the centroid (positions and velocities) 98!> \param pint_env ... 99!> \par History 100!> various bug fixes - hforbert 101!> 2010-11-25 rewritten, added support for velocity printing, 102!> calc of the stddev of the beads turned off [lwalewski] 103!> \author fawzi 104! ************************************************************************************************** 105 SUBROUTINE pint_write_centroids(pint_env) 106 TYPE(pint_env_type), POINTER :: pint_env 107 108 CHARACTER(len=*), PARAMETER :: routineN = 'pint_write_centroids', & 109 routineP = moduleN//':'//routineN 110 INTEGER, PARAMETER :: n_ids = 2, pos_id = 1, vel_id = 2 111 112 CHARACTER(len=default_string_length) :: ext, form, my_middle_name, unit_str 113 CHARACTER(len=default_string_length), DIMENSION(2) :: content_id, middle_name, sect_path, title 114 INTEGER :: handle, handle1, iat, ib, id, idim, & 115 idir, ierr, outformat, should_output, & 116 unit_nr 117 LOGICAL :: new_file 118 REAL(kind=dp) :: nb, ss, unit_conv, vv 119 TYPE(cell_type), POINTER :: cell 120 TYPE(cp_logger_type), POINTER :: logger 121 TYPE(cp_subsys_type), POINTER :: subsys 122 TYPE(f_env_type), POINTER :: f_env 123 TYPE(particle_list_type), POINTER :: particles 124 TYPE(section_vals_type), POINTER :: print_key 125 126 CALL timeset(routineN, handle1) 127 128 CPASSERT(ASSOCIATED(pint_env)) 129 CPASSERT(pint_env%ref_count > 0) 130 131 sect_path(pos_id) = "MOTION%PINT%PRINT%CENTROID_POS" 132 sect_path(vel_id) = "MOTION%PINT%PRINT%CENTROID_VEL" 133 middle_name(pos_id) = "centroid-pos" 134 middle_name(vel_id) = "centroid-vel" 135 content_id(pos_id) = "POS" 136 content_id(vel_id) = "VEL" 137 WRITE (UNIT=title(pos_id), FMT="(A,I8,A,F20.10)") & 138 " i =", pint_env%iter, & 139 ", E =", SUM(pint_env%e_pot_bead)*pint_env%propagator%physpotscale 140 WRITE (UNIT=title(vel_id), FMT="(A,I8,A,F20.10,A,F20.10)") & 141 " i =", pint_env%iter, & 142 ", E_trm =", pint_env%energy(e_kin_thermo_id), & 143 ", E_vir =", pint_env%energy(e_kin_virial_id) 144 145 NULLIFY (logger) 146 logger => cp_get_default_logger() 147 148 CALL pint_u2x(pint_env, ux=pint_env%uv, x=pint_env%v) 149 150 ! iterate over the properties that we know how to print 151 ! (currently positions and velocities) 152 DO id = 1, n_ids 153 154 print_key => section_vals_get_subs_vals(pint_env%input, & 155 TRIM(sect_path(id))) 156 157 should_output = cp_print_key_should_output( & 158 iteration_info=logger%iter_info, & 159 basis_section=print_key) 160 IF (.NOT. BTEST(should_output, cp_p_file)) CONTINUE 161 162 ! get units of measure for output (if available) 163 CALL section_vals_val_get(print_key, "UNIT", & 164 c_val=unit_str) 165 unit_conv = cp_unit_from_cp2k(1.0_dp, TRIM(unit_str)) 166 167 ! get the format for output 168 CALL section_vals_val_get(print_key, "FORMAT", i_val=outformat) 169 170 SELECT CASE (outformat) 171 CASE (dump_dcd, dump_dcd_aligned_cell) 172 form = "UNFORMATTED" 173 ext = ".dcd" 174 CASE (dump_atomic) 175 form = "FORMATTED" 176 ext = "" 177 CASE (dump_xmol) 178 form = "FORMATTED" 179 ext = ".xyz" 180 CASE default 181 CPABORT("") 182 END SELECT 183 184 NULLIFY (f_env, cell, subsys) 185 CALL f_env_add_defaults(f_env_id=pint_env%replicas%f_env_id, & 186 f_env=f_env, handle=handle) 187 CALL force_env_get(force_env=f_env%force_env, & 188 cell=cell, subsys=subsys) 189 CALL cp_subsys_get(subsys, particles=particles) 190 191 ! calculate and copy the requested property 192 ! to the particles structure 193 nb = REAL(pint_env%p, dp) 194 idim = 0 195 DO iat = 1, pint_env%ndim/3 196 DO idir = 1, 3 197 idim = idim + 1 198 ss = 0.0_dp 199 vv = 0.0_dp 200! ss2=0.0_dp 201 DO ib = 1, pint_env%p 202 ss = ss + pint_env%x(ib, idim) 203 vv = vv + pint_env%v(ib, idim) 204! ss2=ss2+pint_env%x(ib,idim)**2 205 END DO 206 particles%els(iat)%r(idir) = ss/nb 207 particles%els(iat)%v(idir) = vv/nb 208! particles%els(iat)%v(idir)=SQRT(ss2/nb-(ss/nb)**2) 209 END DO 210 END DO 211 212 ! set up the output unit number and file name 213 ! for the current property 214 my_middle_name = TRIM(middle_name(id)) 215 unit_nr = cp_print_key_unit_nr(logger=logger, & 216 basis_section=print_key, print_key_path="", & 217 extension=TRIM(ext), middle_name=TRIM(my_middle_name), & 218 local=.FALSE., file_form=form, is_new_file=new_file) 219 220 ! don't write the 0-th frame if the file already exists 221 IF (.NOT. new_file .AND. (pint_env%iter .LE. pint_env%first_step)) THEN 222 CALL cp_print_key_finished_output(unit_nr, logger, & 223 print_key) 224 CONTINUE 225 END IF 226 227 ! actually perform the i/o - on the ionode only 228 IF (unit_nr > 0) THEN 229 230 CALL write_particle_coordinates( & 231 particles%els, & 232 iunit=unit_nr, & 233 output_format=outformat, & 234 content=content_id(id), & 235 title=title(id), & 236 cell=cell, & 237 unit_conv=unit_conv) 238 239 CALL cp_print_key_finished_output(unit_nr, logger, & 240 print_key, "", local=.FALSE.) 241 242 END IF 243 244 CALL f_env_rm_defaults(f_env, ierr, handle) 245 CPASSERT(ierr == 0) 246 247 END DO 248 249 CALL timestop(handle1) 250 RETURN 251 END SUBROUTINE pint_write_centroids 252 253! *************************************************************************** 254!> \brief Write out the trajectory of the beads (positions and velocities) 255!> \param pint_env ... 256!> \par History 257!> 2010-11-25 added support for velocity printing [lwalewski] 258!> \author hforbert 259! ************************************************************************************************** 260 SUBROUTINE pint_write_trajectory(pint_env) 261 TYPE(pint_env_type), POINTER :: pint_env 262 263 CHARACTER(len=*), PARAMETER :: routineN = 'pint_write_trajectory', & 264 routineP = moduleN//':'//routineN 265 INTEGER, PARAMETER :: force_id = 3, n_ids = 3, pos_id = 1, & 266 vel_id = 2 267 268 CHARACTER(len=default_string_length) :: ext, form, ib_str, my_middle_name, & 269 title, unit_str 270 CHARACTER(len=default_string_length), DIMENSION(3) :: content_id, middle_name, sect_path 271 INTEGER :: handle, handle1, iat, ib, id, idim, & 272 idir, ierr, imag_stride, outformat, & 273 should_output, unit_nr 274 LOGICAL :: new_file 275 REAL(kind=dp) :: unit_conv 276 TYPE(cell_type), POINTER :: cell 277 TYPE(cp_logger_type), POINTER :: logger 278 TYPE(cp_subsys_type), POINTER :: subsys 279 TYPE(f_env_type), POINTER :: f_env 280 TYPE(particle_list_type), POINTER :: particles 281 TYPE(section_vals_type), POINTER :: print_key 282 283 CALL timeset(routineN, handle1) 284 285 CPASSERT(ASSOCIATED(pint_env)) 286 CPASSERT(pint_env%ref_count > 0) 287 288 sect_path(pos_id) = "MOTION%PRINT%TRAJECTORY" 289 sect_path(vel_id) = "MOTION%PRINT%VELOCITIES" 290 sect_path(force_id) = "MOTION%PRINT%FORCES" 291 middle_name(pos_id) = "pos-" 292 middle_name(vel_id) = "vel-" 293 middle_name(force_id) = "force-" 294 content_id(pos_id) = "POS" 295 content_id(vel_id) = "VEL" 296 content_id(force_id) = "FORCE" 297 298 NULLIFY (logger) 299 logger => cp_get_default_logger() 300 301 CALL pint_u2x(pint_env, ux=pint_env%uv, x=pint_env%v) 302 303 ! iterate over the properties that we know how to print 304 ! (currently positions and velocities) 305 DO id = 1, n_ids 306 307 print_key => section_vals_get_subs_vals(pint_env%input, & 308 TRIM(sect_path(id))) 309 310 should_output = cp_print_key_should_output( & 311 iteration_info=logger%iter_info, & 312 basis_section=print_key) 313 IF (.NOT. BTEST(should_output, cp_p_file)) CONTINUE 314 315 ! get units of measure for output (if available) 316 CALL section_vals_val_get(print_key, "UNIT", & 317 c_val=unit_str) 318 unit_conv = cp_unit_from_cp2k(1.0_dp, TRIM(unit_str)) 319 320 ! get the format for output 321 CALL section_vals_val_get(print_key, "FORMAT", i_val=outformat) 322 323 SELECT CASE (outformat) 324 CASE (dump_dcd, dump_dcd_aligned_cell) 325 form = "UNFORMATTED" 326 ext = ".dcd" 327 CASE (dump_atomic) 328 form = "FORMATTED" 329 ext = "" 330 CASE (dump_xmol) 331 form = "FORMATTED" 332 ext = ".xyz" 333 CASE default 334 CPABORT("") 335 END SELECT 336 337 NULLIFY (f_env, cell, subsys) 338 CALL f_env_add_defaults(f_env_id=pint_env%replicas%f_env_id, & 339 f_env=f_env, handle=handle) 340 CALL force_env_get(force_env=f_env%force_env, & 341 cell=cell, subsys=subsys) 342 CALL cp_subsys_get(subsys, particles=particles) 343 344 !Get print stride for bead trajectories 345 CALL section_vals_val_get(pint_env%input, & 346 "MOTION%PINT%PRINT%IMAGINARY_TIME_STRIDE", & 347 i_val=imag_stride) 348 349 ! iterate over beads 350 DO ib = 1, pint_env%p, imag_stride 351 352 ! copy the requested property of the current bead 353 ! to the particles structure 354 idim = 0 355 DO iat = 1, pint_env%ndim/3 356 DO idir = 1, 3 357 idim = idim + 1 358 particles%els(iat)%r(idir) = pint_env%x(ib, idim) 359 particles%els(iat)%v(idir) = pint_env%v(ib, idim) 360 particles%els(iat)%f(idir) = pint_env%f(ib, idim) 361 END DO 362 END DO 363 364 ! set up the output unit number and file name 365 ! for the current property and bead 366 ib_str = "" 367 WRITE (ib_str, *) ib 368 my_middle_name = TRIM(middle_name(id))//TRIM(ADJUSTL(ib_str)) 369 unit_nr = cp_print_key_unit_nr(logger=logger, & 370 basis_section=print_key, print_key_path="", & 371 extension=TRIM(ext), middle_name=TRIM(my_middle_name), & 372 local=.FALSE., file_form=form, is_new_file=new_file) 373 374 ! don't write the 0-th frame if the file already exists 375 IF (.NOT. new_file .AND. (pint_env%iter .LE. pint_env%first_step)) THEN 376 CALL cp_print_key_finished_output(unit_nr, logger, & 377 print_key) 378 CONTINUE 379 END IF 380 381 ! actually perform the i/o - on the ionode only 382 IF (unit_nr > 0) THEN 383 384 IF (outformat == dump_xmol) THEN 385 WRITE (UNIT=title, FMT="(A,I8,A,F20.10)") & 386 " i =", pint_env%iter, & 387 ", E =", pint_env%e_pot_bead(ib) 388 END IF 389 390 CALL write_particle_coordinates( & 391 particles%els, & 392 iunit=unit_nr, & 393 output_format=outformat, & 394 content=content_id(id), & 395 title=title, & 396 cell=cell, & 397 unit_conv=unit_conv) 398 399 CALL cp_print_key_finished_output(unit_nr, logger, & 400 print_key, "", local=.FALSE.) 401 402 END IF 403 404 END DO 405 406 CALL f_env_rm_defaults(f_env, ierr, handle) 407 CPASSERT(ierr == 0) 408 409 END DO 410 411 CALL timestop(handle1) 412 RETURN 413 END SUBROUTINE pint_write_trajectory 414 415! *************************************************************************** 416!> \brief Write center of mass (COM) position according to PINT%PRINT%COM 417!> \param pint_env ... 418!> \date 2010-02-17 419!> \author Lukasz Walewski 420! ************************************************************************************************** 421 SUBROUTINE pint_write_com(pint_env) 422 423 TYPE(pint_env_type), POINTER :: pint_env 424 425 CHARACTER(len=*), PARAMETER :: routineN = 'pint_write_com', routineP = moduleN//':'//routineN 426 427 CHARACTER(len=default_string_length) :: stmp1, stmp2 428 INTEGER :: ic, unit_nr 429 LOGICAL :: new_file, should_output 430 REAL(kind=dp), DIMENSION(3) :: com_r 431 TYPE(cp_logger_type), POINTER :: logger 432 TYPE(section_vals_type), POINTER :: print_key 433 434 CPASSERT(ASSOCIATED(pint_env)) 435 436 NULLIFY (logger) 437 logger => cp_get_default_logger() 438 439 ! decide whether to write anything or not 440 NULLIFY (print_key) 441 print_key => section_vals_get_subs_vals(pint_env%input, & 442 "MOTION%PINT%PRINT%COM") 443 should_output = BTEST(cp_print_key_should_output( & 444 iteration_info=logger%iter_info, & 445 basis_section=print_key), cp_p_file) 446 IF (.NOT. should_output) THEN 447 RETURN 448 END IF 449 450 com_r = pint_com_pos(pint_env) 451 DO ic = 1, 3 452 com_r(ic) = cp_unit_from_cp2k(com_r(ic), "angstrom") 453 END DO 454 455 unit_nr = cp_print_key_unit_nr(logger, print_key, is_new_file=new_file, & 456 middle_name="com-pos", extension=".xyz") 457 458 ! don't write the 0-th frame if the file already exists 459 IF (.NOT. new_file .AND. (pint_env%iter .LE. pint_env%first_step)) THEN 460 CALL cp_print_key_finished_output(unit_nr, logger, & 461 print_key) 462 RETURN 463 END IF 464 465 ! actually perform the i/o - on the ionode only 466 IF (unit_nr > 0) THEN 467 468 WRITE (unit_nr, '(I2)') 1 469 WRITE (stmp1, *) pint_env%iter 470 WRITE (stmp2, '(F20.10)') pint_env%energy(e_conserved_id) 471 WRITE (unit_nr, '(4A)') " Iteration = ", TRIM(ADJUSTL(stmp1)), & 472 ", E_conserved = ", TRIM(ADJUSTL(stmp2)) 473 WRITE (unit_nr, '(A2,3(1X,F20.10))') "X ", (com_r(ic), ic=1, 3) 474 475 CALL m_flush(unit_nr) 476 477 END IF 478 479 CALL cp_print_key_finished_output(unit_nr, logger, print_key) 480 481 RETURN 482 END SUBROUTINE pint_write_com 483 484! *************************************************************************** 485!> \brief Writes out the energies according to PINT%PRINT%ENERGY 486!> \param pint_env path integral environment 487!> \par History 488!> various bug fixes [hforbert] 489!> 2009-11-16 energy components calc moved out of here [lwalewski] 490!> \author fawzi 491! ************************************************************************************************** 492 SUBROUTINE pint_write_ener(pint_env) 493 TYPE(pint_env_type), POINTER :: pint_env 494 495 CHARACTER(len=*), PARAMETER :: routineN = 'pint_write_ener', & 496 routineP = moduleN//':'//routineN 497 498 INTEGER :: ndof, unit_nr 499 LOGICAL :: file_is_new 500 REAL(kind=dp) :: t, temp 501 TYPE(cp_logger_type), POINTER :: logger 502 TYPE(section_vals_type), POINTER :: print_key 503 504 CPASSERT(ASSOCIATED(pint_env)) 505 CPASSERT(pint_env%ref_count > 0) 506 507 NULLIFY (print_key, logger) 508 print_key => section_vals_get_subs_vals(pint_env%input, & 509 "MOTION%PINT%PRINT%ENERGY") 510 logger => cp_get_default_logger() 511 IF (BTEST(cp_print_key_should_output(iteration_info=logger%iter_info, & 512 basis_section=print_key), cp_p_file)) THEN 513 514 unit_nr = cp_print_key_unit_nr(logger, print_key, middle_name="energy", & 515 extension=".dat", is_new_file=file_is_new) 516 517 ! don't write the 0-th frame if the file already exists 518 IF (.NOT. file_is_new .AND. (pint_env%iter .LE. pint_env%first_step)) THEN 519 CALL cp_print_key_finished_output(unit_nr, logger, & 520 print_key) 521 RETURN 522 END IF 523 524 ! cp_print_key_unit_nr returns -1 on nodes other than logger%para_env%ionode 525 IF (unit_nr > 0) THEN 526 527 ! please keep the format explanation up to date 528 ! keep the constant of motion the true constant of motion ! 529 IF (file_is_new) THEN 530 WRITE (unit_nr, "(A8,1X,A12,1X,5(A20,1X),A12)") & 531 "# StepNr", & 532 " Time [fs]", & 533 " Kinetic [a.u.]", & 534 " VirialKin [a.u.]", & 535 " Temperature [K]", & 536 " Potential [a.u.]", & 537 " ConsQty [a.u.]", & 538 " CPU [s]" 539 END IF 540 541 t = cp_unit_from_cp2k(pint_env%t, "fs") 542 543 ndof = pint_env%p 544 IF (pint_env%first_propagated_mode .EQ. 2) THEN 545 ndof = ndof - 1 546 END IF 547 temp = cp_unit_from_cp2k(2.0_dp*pint_env%e_kin_beads/ & 548 REAL(ndof, dp)/REAL(pint_env%ndim, dp), & 549 "K")*pint_env%propagator%temp_sim2phys 550 551 WRITE (unit_nr, "(I8,1X,F12.3,1X,5(F20.9,1X),F12.1)") & 552 pint_env%iter, & 553 t, & 554 pint_env%energy(e_kin_thermo_id), & 555 pint_env%energy(e_kin_virial_id), & 556 temp, & 557 pint_env%energy(e_potential_id), & 558 pint_env%energy(e_conserved_id), & 559 pint_env%time_per_step 560 CALL m_flush(unit_nr) 561 562 END IF 563 564 CALL cp_print_key_finished_output(unit_nr, logger, print_key) 565 END IF 566 567 RETURN 568 END SUBROUTINE pint_write_ener 569 570! *************************************************************************** 571!> \brief Writes out the actions according to PINT%PRINT%ACTION 572!> \param pint_env path integral environment 573!> \author Felix Uhl 574! ************************************************************************************************** 575 SUBROUTINE pint_write_action(pint_env) 576 TYPE(pint_env_type), POINTER :: pint_env 577 578 CHARACTER(len=*), PARAMETER :: routineN = 'pint_write_action', & 579 routineP = moduleN//':'//routineN 580 581 INTEGER :: unit_nr 582 LOGICAL :: file_is_new 583 REAL(kind=dp) :: t 584 TYPE(cp_logger_type), POINTER :: logger 585 TYPE(section_vals_type), POINTER :: print_key 586 587 CPASSERT(ASSOCIATED(pint_env)) 588 CPASSERT(pint_env%ref_count > 0) 589 590 NULLIFY (print_key, logger) 591 print_key => section_vals_get_subs_vals(pint_env%input, & 592 "MOTION%PINT%PRINT%ACTION") 593 logger => cp_get_default_logger() 594 IF (BTEST(cp_print_key_should_output(iteration_info=logger%iter_info, & 595 basis_section=print_key), cp_p_file)) THEN 596 597 unit_nr = cp_print_key_unit_nr(logger, print_key, middle_name="action", & 598 extension=".dat", is_new_file=file_is_new) 599 600 ! don't write the 0-th frame if the file already exists 601 IF (.NOT. file_is_new .AND. (pint_env%iter .LE. pint_env%first_step)) THEN 602 CALL cp_print_key_finished_output(unit_nr, logger, & 603 print_key) 604 RETURN 605 END IF 606 607 ! cp_print_key_unit_nr returns -1 on nodes other than logger%para_env%ionode 608 IF (unit_nr > 0) THEN 609 610 ! please keep the format explanation up to date 611 ! keep the constant of motion the true constant of motion ! 612 IF (file_is_new) THEN 613 WRITE (unit_nr, "(A8,1X,A12,1X,2(A25,1X))") & 614 "# StepNr", & 615 " Time [fs]", & 616 " Link Action [a.u.]", & 617 " Potential Action [a.u.]" 618 END IF 619 620 t = cp_unit_from_cp2k(pint_env%t, "fs") 621 622 WRITE (unit_nr, "(I8,1X,F12.3,1X,5(F20.9,1X),F12.1)") & 623 pint_env%iter, & 624 t, & 625 pint_env%link_action, & 626 pint_env%pot_action 627 CALL m_flush(unit_nr) 628 629 END IF 630 631 CALL cp_print_key_finished_output(unit_nr, logger, print_key) 632 END IF 633 634 RETURN 635 END SUBROUTINE pint_write_action 636 637 ! *************************************************************************** 638 !> \brief Write step info to the output file. 639 !> \param pint_env ... 640 !> \date 2009-11-16 641 !> \par History 642 !> 2010-01-27 getting default unit nr now only on ionode [lwalewski] 643 !> \author Lukasz Walewski 644 ! ************************************************************************************************** 645! ************************************************************************************************** 646!> \brief ... 647!> \param pint_env ... 648! ************************************************************************************************** 649 SUBROUTINE pint_write_step_info(pint_env) 650 TYPE(pint_env_type), POINTER :: pint_env 651 652 CHARACTER(len=*), PARAMETER :: routineN = 'pint_write_step_info', & 653 routineP = moduleN//':'//routineN 654 655 CHARACTER(len=default_string_length) :: msgstr, stmp, time_unit 656 INTEGER :: unit_nr 657 REAL(kind=dp) :: time_used 658 TYPE(cp_logger_type), POINTER :: logger 659 660 CPASSERT(ASSOCIATED(pint_env)) 661 662 unit_nr = 0 663 NULLIFY (logger) 664 logger => cp_get_default_logger() 665 666 time_used = pint_env%time_per_step 667 time_unit = "sec" 668 IF (time_used .GE. 60.0_dp) THEN 669 time_used = time_used/60.0_dp 670 time_unit = "min" 671 END IF 672 IF (time_used .GE. 60.0_dp) THEN 673 time_used = time_used/60.0_dp 674 time_unit = "hours" 675 END IF 676 msgstr = "PINT step" 677 stmp = "" 678 WRITE (stmp, *) pint_env%iter 679 msgstr = TRIM(ADJUSTL(msgstr))//" "//TRIM(ADJUSTL(stmp))//" of" 680 stmp = "" 681 WRITE (stmp, *) pint_env%last_step 682 msgstr = TRIM(ADJUSTL(msgstr))//" "//TRIM(ADJUSTL(stmp))//" in" 683 stmp = "" 684 WRITE (stmp, '(F20.1)') time_used 685 msgstr = TRIM(ADJUSTL(msgstr))//" "//TRIM(ADJUSTL(stmp)) 686 msgstr = TRIM(ADJUSTL(msgstr))//" "//TRIM(ADJUSTL(time_unit))//"." 687 688 IF (logger%para_env%ionode) THEN 689 unit_nr = cp_logger_get_default_unit_nr(logger) 690 WRITE (unit_nr, '(T2,A)') "PINT| "//TRIM(ADJUSTL(msgstr)) 691 END IF 692 693 ! print out the total energy - for regtest evaluation 694 stmp = "" 695 WRITE (stmp, *) pint_env%energy(e_conserved_id) 696 msgstr = "Total energy = "//TRIM(ADJUSTL(stmp)) 697 IF (logger%para_env%ionode) THEN 698 WRITE (unit_nr, '(T2,A)') "PINT| "//TRIM(ADJUSTL(msgstr)) 699 END IF 700 701 RETURN 702 END SUBROUTINE pint_write_step_info 703 704! *************************************************************************** 705!> \brief Write radii of gyration according to PINT%PRINT%CENTROID_GYR 706!> \param pint_env ... 707!> \date 2011-01-07 708!> \author Lukasz Walewski 709! ************************************************************************************************** 710 SUBROUTINE pint_write_rgyr(pint_env) 711 712 TYPE(pint_env_type), POINTER :: pint_env 713 714 CHARACTER(len=*), PARAMETER :: routineN = 'pint_write_rgyr', & 715 routineP = moduleN//':'//routineN 716 717 CHARACTER(len=default_string_length) :: unit_str 718 INTEGER :: ia, ib, ic, idim, unit_nr 719 LOGICAL :: new_file, should_output 720 REAL(kind=dp) :: nb, ss, unit_conv 721 TYPE(cp_logger_type), POINTER :: logger 722 TYPE(section_vals_type), POINTER :: print_key 723 724 CPASSERT(ASSOCIATED(pint_env)) 725 726 NULLIFY (logger) 727 logger => cp_get_default_logger() 728 729 ! decide whether to write anything or not 730 NULLIFY (print_key) 731 print_key => section_vals_get_subs_vals(pint_env%input, & 732 "MOTION%PINT%PRINT%CENTROID_GYR") 733 should_output = BTEST(cp_print_key_should_output( & 734 iteration_info=logger%iter_info, & 735 basis_section=print_key), cp_p_file) 736 IF (.NOT. should_output) THEN 737 RETURN 738 END IF 739 740 ! get the units conversion factor 741 CALL section_vals_val_get(print_key, "UNIT", c_val=unit_str) 742 unit_conv = cp_unit_from_cp2k(1.0_dp, TRIM(unit_str)) 743 744 ! calculate the centroid positions 745 nb = REAL(pint_env%p, dp) 746 idim = 0 747 DO ia = 1, pint_env%ndim/3 748 DO ic = 1, 3 749 idim = idim + 1 750 ss = 0.0_dp 751 DO ib = 1, pint_env%p 752 ss = ss + pint_env%x(ib, idim) 753 END DO 754 pint_env%rtmp_ndim(idim) = ss/nb 755 END DO 756 END DO 757 758 ! calculate the radii of gyration 759 idim = 0 760 DO ia = 1, pint_env%ndim/3 761 ss = 0.0_dp 762 DO ic = 1, 3 763 idim = idim + 1 764 DO ib = 1, pint_env%p 765 ss = ss + (pint_env%x(ib, idim) - pint_env%rtmp_ndim(idim))**2 766 END DO 767 END DO 768 pint_env%rtmp_natom(ia) = SQRT(ss/nb)*unit_conv 769 END DO 770 771 unit_nr = cp_print_key_unit_nr(logger, print_key, is_new_file=new_file, & 772 middle_name="centroid-gyr", extension=".dat") 773 774 ! don't write the 0-th frame if the file already exists 775 IF (.NOT. new_file .AND. (pint_env%iter .LE. pint_env%first_step)) THEN 776 CALL cp_print_key_finished_output(unit_nr, logger, & 777 print_key) 778 RETURN 779 END IF 780 781 ! actually perform the i/o - on the ionode only 782 IF (unit_nr > 0) THEN 783 784 DO ia = 1, pint_env%ndim/3 785 WRITE (unit_nr, '(F20.10,1X)', ADVANCE='NO') pint_env%rtmp_natom(ia) 786 END DO 787 WRITE (unit_nr, '(A)') "" 788 789 CALL m_flush(unit_nr) 790 791 END IF 792 793 CALL cp_print_key_finished_output(unit_nr, logger, print_key) 794 795 RETURN 796 END SUBROUTINE pint_write_rgyr 797 798END MODULE pint_io 799