TABLE OF CONTENTS
- ABINIT/m_rttddft_output
- m_rttddft_output/prt_den
- m_rttddft_output/prt_dos
- m_rttddft_output/prt_eig
- m_rttddft_output/prt_occ
- m_rttddft_output/prt_restart
- m_rttddft_output/prt_wfk
- m_rttddft_output/rttddft_output
ABINIT/m_rttddft_output [ Modules ]
NAME
m_rttddft_ouptut
FUNCTION
Manages most output of RT-TDDFT runs
COPYRIGHT
Copyright (C) 2021-2024 ABINIT group (FB) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_rttddft_output 23 24 use defs_basis 25 use defs_abitypes, only: MPI_type 26 use defs_datatypes, only: pseudopotential_type, ebands_t 27 28 use m_common, only: prteigrs 29 use m_crystal, only: crystal_init, crystal_t 30 use m_dtfil, only: datafiles_type 31 use m_dtset, only: dataset_type 32 use m_ebands, only: ebands_init, ebands_free 33 use m_epjdos, only: dos_calcnwrite, partial_dos_fractions, & 34 & partial_dos_fractions_paw,epjdos_t, & 35 & epjdos_new, prtfatbands, fatbands_ncwrite 36 use m_errors, only: msg_hndl, assert 37 use m_ioarr, only: fftdatar_write 38 use m_io_tools, only: open_file 39 use m_iowf, only: outwf 40 use m_mpinfo, only: iwrite_fftdatar 41 #ifdef HAVE_NETCDF 42 use netcdf 43 #endif 44 use m_paral_atom, only: get_my_atmtab, free_my_atmtab 45 use m_profiling_abi, only: abimem_record 46 use m_rttddft_tdks, only: tdks_type 47 use m_specialmsg, only: wrtout 48 use m_xmpi, only: xmpi_comm_rank 49 50 implicit none 51 52 private
m_rttddft_output/prt_den [ Functions ]
[ Top ] [ m_rttddft_output ] [ Functions ]
NAME
prt_den
FUNCTION
Outputs the electronic density
INPUTS
dtfil <type datafiles_type> = infos about file names, file unit numbers dtset <type(dataset_type)> = all input variables for this dataset istep <integer> = step number mpi_enreg <MPI_type> = MPI-parallelisation information psps <type(pseudopotential_type)> = variables related to pseudopotentials tdks <type(tdks_type)> = the tdks object to initialize
OUTPUT
SOURCE
372 subroutine prt_den(dtfil, dtset, istep, mpi_enreg, psps, tdks) 373 374 implicit none 375 376 !Arguments ------------------------------------ 377 !scalars 378 integer, intent(in) :: istep 379 type(datafiles_type), intent(inout) :: dtfil 380 type(dataset_type), intent(inout) :: dtset 381 type(MPI_type), intent(inout) :: mpi_enreg 382 type(pseudopotential_type), intent(inout) :: psps 383 type(tdks_type), intent(inout) :: tdks 384 !arrays 385 386 !Local variables------------------------------- 387 !scalars 388 integer,parameter :: cplex1=1 389 integer :: bantot 390 integer :: iatom 391 integer :: spacecomm 392 integer :: my_comm_atom, my_natom 393 integer :: me 394 integer :: natom 395 !#ifdef HAVE_NETCDF 396 ! integer :: ncid 397 !#endif 398 integer :: timrev 399 character(len=fnlen) :: fname 400 character(len=24) :: step_nb 401 logical :: paral_atom 402 logical :: remove_inv 403 logical :: my_atmtab_allocated 404 type(crystal_t) :: crystal 405 type(ebands_t) :: ebands 406 !arrays 407 integer, pointer :: my_atmtab(:) 408 real(dp), allocatable :: doccde(:) 409 410 ! ************************************************************************* 411 412 if (dtset%prtden /= 0) then 413 spacecomm = mpi_enreg%comm_cell 414 me = xmpi_comm_rank(spacecomm) 415 416 natom = dtset%natom 417 my_natom = mpi_enreg%my_natom 418 paral_atom=(my_natom/=natom) 419 my_comm_atom = mpi_enreg%comm_atom 420 nullify(my_atmtab) 421 if (paral_atom) then 422 call get_my_atmtab(mpi_enreg%comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 423 else 424 ABI_MALLOC(my_atmtab, (natom)) 425 my_atmtab = (/ (iatom, iatom=1, natom) /) 426 my_atmtab_allocated = .true. 427 end if 428 429 !FB: Maybe this should be moved out of that subroutine if needed in other outputs than densities 430 remove_inv=.false. 431 timrev = 2; if (any(dtset%kptopt == [3, 4])) timrev= 1 432 call crystal_init(dtset%amu_orig(:,1),crystal,dtset%spgroup,natom,dtset%npsp,psps%ntypat, & 433 dtset%nsym,tdks%rprimd,dtset%typat,tdks%xred,dtset%ziontypat,dtset%znucl,timrev,& 434 dtset%nspden==2.and.dtset%nsppol==1,remove_inv,tdks%hdr%title,& 435 dtset%symrel,dtset%tnons,dtset%symafm) 436 !Electron band energies. 437 bantot= dtset%mband*dtset%nkpt*dtset%nsppol 438 ABI_CALLOC(doccde, (bantot)) 439 call ebands_init(bantot,ebands,dtset%nelect,dtset%ne_qFD,dtset%nh_qFD,dtset%ivalence, & 440 doccde,tdks%eigen,dtset%istwfk,dtset%kptns,dtset%nband,dtset%nkpt,tdks%npwarr,dtset%nsppol, & 441 dtset%nspinor,dtset%tphysel,dtset%tsmear,dtset%occopt,tdks%occ0,dtset%wtk,& 442 dtset%cellcharge(1),dtset%kptopt,dtset%kptrlatt_orig,dtset%nshiftk_orig,dtset%shiftk_orig, & 443 dtset%kptrlatt,dtset%nshiftk,dtset%shiftk) 444 ABI_FREE(doccde) 445 446 write(step_nb,*) istep 447 448 !** Outputs the density 449 !Warnings : 450 !- core charge is excluded from the charge density; 451 !- the potential is the INPUT vtrial. 452 if (iwrite_fftdatar(mpi_enreg)) then 453 fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_DEN' 454 call fftdatar_write("density",fname,dtset%iomode,tdks%hdr,crystal,tdks%pawfgr%ngfft, & 455 & cplex1,tdks%pawfgr%nfft,dtset%nspden,tdks%rhor,mpi_enreg,ebands=ebands) 456 end if 457 458 call crystal%free() 459 call ebands_free(ebands) 460 end if 461 462 end subroutine prt_den
m_rttddft_output/prt_dos [ Functions ]
[ Top ] [ m_rttddft_output ] [ Functions ]
NAME
prt_dos
FUNCTION
Computes and outputs the electronic DOS
INPUTS
dtfil <type datafiles_type> = infos about file names, file unit numbers dtset <type(dataset_type)> = all input variables for this dataset istep <integer> = step number mpi_enreg <MPI_type> = MPI-parallelisation information psps <type(pseudopotential_type)> = variables related to pseudopotentials tdks <type(tdks_type)> = the tdks object to initialize
OUTPUT
SOURCE
484 subroutine prt_dos(dtfil, dtset, istep, mpi_enreg, psps, tdks) 485 486 implicit none 487 488 !Arguments ------------------------------------ 489 !scalars 490 integer, intent(in) :: istep 491 type(datafiles_type), intent(inout) :: dtfil 492 type(dataset_type), intent(inout) :: dtset 493 type(MPI_type), intent(inout) :: mpi_enreg 494 type(pseudopotential_type), intent(inout) :: psps 495 type(tdks_type), intent(inout) :: tdks 496 !arrays 497 498 !Local variables------------------------------- 499 !scalars 500 integer,parameter :: master=0 501 integer :: bantot 502 integer :: collect 503 integer :: iatom 504 integer :: spacecomm 505 integer :: my_comm_atom, my_natom 506 integer :: me 507 integer :: natom 508 !#ifdef HAVE_NETCDF 509 ! integer :: ncid 510 !#endif 511 integer :: timrev 512 character(len=fnlen) :: fname 513 character(len=24) :: step_nb 514 logical :: paral_atom 515 logical :: remove_inv 516 logical :: my_atmtab_allocated 517 type(crystal_t) :: crystal 518 type(epjdos_t) :: dos 519 type(ebands_t) :: ebands 520 !arrays 521 integer, pointer :: my_atmtab(:) 522 real(dp), allocatable :: doccde(:) 523 524 ! ************************************************************************* 525 526 spacecomm = mpi_enreg%comm_cell 527 me = xmpi_comm_rank(spacecomm) 528 529 !FB: @MT - Is this needed? 530 natom = dtset%natom 531 my_natom = mpi_enreg%my_natom 532 paral_atom=(my_natom/=natom) 533 my_comm_atom = mpi_enreg%comm_atom 534 nullify(my_atmtab) 535 if (paral_atom) then 536 call get_my_atmtab(mpi_enreg%comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 537 else 538 ABI_MALLOC(my_atmtab, (natom)) 539 my_atmtab = (/ (iatom, iatom=1, natom) /) 540 my_atmtab_allocated = .true. 541 end if 542 543 remove_inv=.false. 544 timrev = 2; if (any(dtset%kptopt == [3, 4])) timrev= 1 545 call crystal_init(dtset%amu_orig(:,1),crystal,dtset%spgroup,natom,dtset%npsp,psps%ntypat, & 546 dtset%nsym,tdks%rprimd,dtset%typat,tdks%xred,dtset%ziontypat,dtset%znucl,timrev,& 547 dtset%nspden==2.and.dtset%nsppol==1,remove_inv,tdks%hdr%title,& 548 dtset%symrel,dtset%tnons,dtset%symafm) 549 !Electron band energies. 550 bantot= dtset%mband*dtset%nkpt*dtset%nsppol 551 ABI_CALLOC(doccde, (bantot)) 552 call ebands_init(bantot,ebands,dtset%nelect,dtset%ne_qFD,dtset%nh_qFD,dtset%ivalence, & 553 doccde,tdks%eigen,dtset%istwfk,dtset%kptns,dtset%nband,dtset%nkpt,tdks%npwarr,dtset%nsppol, & 554 dtset%nspinor,dtset%tphysel,dtset%tsmear,dtset%occopt,tdks%occ0,dtset%wtk,& 555 dtset%cellcharge(1),dtset%kptopt,dtset%kptrlatt_orig,dtset%nshiftk_orig,dtset%shiftk_orig, & 556 dtset%kptrlatt,dtset%nshiftk,dtset%shiftk) 557 ABI_FREE(doccde) 558 559 write(step_nb,*) istep 560 561 !** Generate DOS using the tetrahedron method or using Gaussians 562 if (dtset%prtdos>=2.or.dtset%pawfatbnd>0) then 563 dos = epjdos_new(dtset, psps, tdks%pawtab) 564 565 if (dos%partial_dos_flag>=1 .or. dos%fatbands_flag==1)then 566 ! Generate fractions for partial DOSs if needed partial_dos 1,2,3,4 give different decompositions 567 collect = 1 !; if (psps%usepaw==1 .and. dos%partial_dos_flag /= 2) collect = 0 568 if ((psps%usepaw==0.or.dtset%pawprtdos/=2) .and. dos%partial_dos_flag>=1) then 569 call partial_dos_fractions(dos,crystal,dtset,tdks%eigen,tdks%occ0,tdks%npwarr,tdks%kg,tdks%cg,tdks%mcg,collect,mpi_enreg) 570 end if 571 572 if (psps%usepaw==1 .and. dos%partial_dos_flag /= 2) then 573 ! TODO: update partial_dos_fractions_paw for extra atoms - no PAW contribution normally, but check bounds and so on. 574 call partial_dos_fractions_paw(dos,tdks%cprj,tdks%dimcprj,dtset,tdks%mcprj,dtset%mkmem,mpi_enreg,tdks%pawrad,tdks%pawtab) 575 end if 576 else 577 dos%fractions(:,:,:,1)=one 578 end if 579 580 !Here, print out fatbands for the k-points given in file appended _FATBANDS 581 if (me == master .and. dtset%pawfatbnd>0 .and. dos%fatbands_flag==1) then 582 fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_FATBANDS' 583 call prtfatbands(dos,dtset,ebands,fname,dtset%pawfatbnd,tdks%pawtab) 584 end if 585 586 !Here, computation and output of DOS and partial DOS _DOS 587 if (dos%fatbands_flag == 0 .and. dos%prtdos /= 4) then 588 fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_DOS' 589 call dos_calcnwrite(dos,dtset,crystal,ebands,fname,spacecomm) 590 end if 591 end if 592 593 call dos%free() 594 call crystal%free() 595 call ebands_free(ebands) 596 597 end subroutine prt_dos
m_rttddft_output/prt_eig [ Functions ]
[ Top ] [ m_rttddft_output ] [ Functions ]
NAME
prt_eig
FUNCTION
Outputs eigenvalues
INPUTS
dtfil <type datafiles_type> = infos about file names, file unit numbers dtset <type(dataset_type)> = all input variables for this dataset istep <integer> = step number mpi_enreg <MPI_type> = MPI-parallelisation information tdks <type(tdks_type)> = the tdks object to initialize
OUTPUT
SIDE EFFECTS
SOURCE
203 subroutine prt_eig(dtfil, dtset, istep, mpi_enreg, tdks) 204 205 implicit none 206 207 !Arguments ------------------------------------ 208 !scalars 209 integer, intent(in) :: istep 210 type(datafiles_type), intent(inout) :: dtfil 211 type(dataset_type), intent(inout) :: dtset 212 type(MPI_type), intent(inout) :: mpi_enreg 213 type(tdks_type), intent(inout) :: tdks 214 !arrays 215 216 !Local variables------------------------------- 217 !scalars 218 integer,parameter :: enunit=0, option=3 219 integer :: me 220 integer :: spacecomm 221 real(dp) :: vxcavg_dum 222 character(len=fnlen) :: fname 223 character(len=24) :: step_nb 224 !arrays 225 real(dp) :: resid(dtset%mband*dtset%nkpt*dtset%nsppol) 226 227 ! ************************************************************************* 228 229 spacecomm = mpi_enreg%comm_cell 230 me = xmpi_comm_rank(spacecomm) 231 232 write(step_nb,*) istep 233 fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_EIG' 234 resid = zero 235 vxcavg_dum=zero 236 237 if(me==0)then 238 call prteigrs(tdks%eigen,enunit,tdks%energies%e_fermie,tdks%energies%e_fermih, & 239 & fname,ab_out,dtset%iscf,dtset%kptns,dtset%kptopt,dtset%mband, & 240 & dtset%nband,dtset%nbdbuf,dtset%nkpt,0,dtset%nsppol,tdks%occ0, & 241 & dtset%occopt,option,dtset%prteig,dtset%prtvol,resid,dtset%tolwfr, & 242 & vxcavg_dum,dtset%wtk) 243 end if 244 245 end subroutine prt_eig
m_rttddft_output/prt_occ [ Functions ]
[ Top ] [ m_rttddft_output ] [ Functions ]
NAME
prt_occ
FUNCTION
Outputs occupation numbers
INPUTS
dtfil <type datafiles_type> = infos about file names, file unit numbers dtset <type(dataset_type)> = all input variables for this dataset istep <integer> = step number mpi_enreg <MPI_type> = MPI-parallelisation information tdks <type(tdks_type)> = the tdks object to initialize
OUTPUT
SOURCE
266 subroutine prt_occ(dtfil, dtset, istep, mpi_enreg, tdks) 267 268 implicit none 269 270 !Arguments ------------------------------------ 271 !scalars 272 integer, intent(in) :: istep 273 type(datafiles_type), intent(inout) :: dtfil 274 type(dataset_type), intent(inout) :: dtset 275 type(MPI_type), intent(inout) :: mpi_enreg 276 type(tdks_type), intent(inout) :: tdks 277 !arrays 278 279 !Local variables------------------------------- 280 !scalars 281 integer :: band_index 282 integer :: iband, ii, ikpt, isppol 283 integer :: me 284 integer :: nkpt 285 integer :: nband_k, nsppol 286 integer :: temp_unit 287 !arrays 288 character(len=fnlen) :: fname 289 character(len=4) :: ibnd_fmt, ikpt_fmt 290 character(len=500) :: msg 291 character(len=24) :: step_nb 292 293 294 ! ************************************************************************* 295 296 if (dtset%prtocc > 0) then 297 me = xmpi_comm_rank(mpi_enreg%comm_cell) 298 299 write(step_nb,*) istep 300 fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_OCC' 301 302 if (open_file(fname, msg, newunit=temp_unit, status='unknown', form='formatted') /= 0) then 303 ABI_ERROR(msg) 304 end if 305 306 nkpt = dtset%nkpt 307 nsppol = dtset%nsppol 308 309 if(me==0)then 310 band_index=0 311 do isppol=1,nsppol 312 313 if(nsppol==2)then 314 if(isppol==1)write(msg, '(2a)' ) ch10,' SPIN UP channel ' 315 if(isppol==2)write(msg, '(2a)' ) ch10,' SPIN DOWN channel ' 316 call wrtout(temp_unit,msg) 317 end if 318 ikpt_fmt="i4" ; if(nkpt>=10000)ikpt_fmt="i6" ; if(nkpt>=1000000)ikpt_fmt="i9" 319 if (nsppol==2.and.isppol==1) then 320 write(msg, '(a,'//ikpt_fmt//',2x,a)' ) & 321 'Occupation numbers for nkpt=',nkpt,'k points, SPIN UP:' 322 else if (nsppol==2.and.isppol==2) then 323 write(msg, '(a,'//ikpt_fmt//',2x,a)' ) & 324 'Occupation numbers for nkpt=',nkpt,'k points, SPIN DOWN:' 325 else 326 write(msg, '(a,'//ikpt_fmt//',2x,a)' ) & 327 'Occupation numbers for nkpt=',nkpt,'k points:' 328 end if 329 call wrtout(temp_unit,msg) 330 do ikpt=1,nkpt 331 nband_k=dtset%nband(ikpt+(isppol-1)*nkpt) 332 ikpt_fmt="i4" ; if(nkpt>=10000)ikpt_fmt="i6" ; if(nkpt>=1000000)ikpt_fmt="i9" 333 ibnd_fmt="i3" ; if(nband_k>=1000)ibnd_fmt="i6" ; if(nband_k>=1000000)ibnd_fmt="i9" 334 write(msg, '(a,'//ikpt_fmt//',a,'//ibnd_fmt//',a,f9.5,a,3f8.4,a)' ) & 335 & ' kpt#',ikpt,', nband=',nband_k,', wtk=',dtset%wtk(ikpt)+tol10,', kpt=',& 336 & dtset%kptns(1:3,ikpt)+tol10,' (reduced coord)' 337 call wrtout(temp_unit,msg) 338 do ii=0,(nband_k-1)/6 339 write(msg, '(1p,6e12.4)')(tdks%occ(iband+band_index),iband=1+6*ii,min(6+6*ii,nband_k)) 340 call wrtout(temp_unit,msg) 341 end do 342 band_index=band_index+nband_k 343 end do 344 end do 345 end if 346 347 close(temp_unit) 348 end if 349 350 end subroutine prt_occ
m_rttddft_output/prt_restart [ Functions ]
[ Top ] [ m_rttddft_output ] [ Functions ]
NAME
prt_restart
FUNCTION
Print restart file
INPUTS
dtfil <type datafiles_type> = infos about file names, file unit numbers istep <integer> = step number mpi_enreg <MPI_type> = MPI-parallelisation information tdks <type(tdks_type)> = the tdks object to initialize
OUTPUT
SOURCE
685 subroutine prt_restart(dtfil, istep, mpi_enreg, tdks) 686 687 implicit none 688 689 !Arguments ------------------------------------ 690 !scalars 691 integer, intent(in) :: istep 692 type(datafiles_type), intent(inout) :: dtfil 693 type(MPI_type), intent(inout) :: mpi_enreg 694 type(tdks_type), intent(inout) :: tdks 695 !arrays 696 697 !Local variables------------------------------- 698 !scalars 699 character(len=500) :: msg 700 character(len=fnlen) :: fname 701 character(len=24 ) :: step_nb 702 !arrays 703 704 ! ************************************************************************* 705 706 if (mpi_enreg%me == 0) then 707 write(step_nb,*) istep 708 rewind(tdks%tdrestart_unit) 709 write(msg,'(a)') step_nb 710 call wrtout(tdks%tdrestart_unit,msg) 711 write(msg,'(a)') tdks%fname_tdener 712 call wrtout(tdks%tdrestart_unit,msg) 713 write(msg,'(a)') tdks%fname_wfk0 714 call wrtout(tdks%tdrestart_unit,msg) 715 fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_WFK' 716 write(msg,'(a)') fname 717 call wrtout(tdks%tdrestart_unit,msg) 718 end if 719 720 end subroutine prt_restart
m_rttddft_output/prt_wfk [ Functions ]
[ Top ] [ m_rttddft_output ] [ Functions ]
NAME
prt_wfk
FUNCTION
Outputs wavefunctions in WFK file
INPUTS
dtfil <type datafiles_type> = infos about file names, file unit numbers dtset <type(dataset_type)> = all input variables for this dataset istep <integer> = step number mpi_enreg <MPI_type> = MPI-parallelisation information psps <type(pseudopotential_type)> = variables related to pseudopotentials tdks <type(tdks_type)> = the tdks object to initialize force_write <logical> = force the writing of WFK (useful for last step) - optional
OUTPUT
SOURCE
620 subroutine prt_wfk(dtfil, dtset, istep, mpi_enreg, psps, tdks, force_write) 621 622 implicit none 623 624 !Arguments ------------------------------------ 625 !scalars 626 integer, intent(in) :: istep 627 type(datafiles_type), intent(inout) :: dtfil 628 type(dataset_type), intent(inout) :: dtset 629 type(MPI_type), intent(inout) :: mpi_enreg 630 type(pseudopotential_type), intent(inout) :: psps 631 type(tdks_type), intent(inout) :: tdks 632 logical, optional, intent(in) :: force_write 633 !arrays 634 635 !Local variables------------------------------- 636 !scalars 637 integer,parameter :: response=0 638 character(len=fnlen) :: fname 639 character(len=24) :: step_nb 640 logical :: lforce_write = .FALSE. 641 !arrays 642 643 ! ************************************************************************* 644 645 write(step_nb,*) istep 646 fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_WFK' 647 648 if (present(force_write)) then 649 if (force_write) lforce_write = .TRUE. 650 end if 651 652 !Use initial eigenvalues to ensure that we get the same occupation upon restart 653 if (lforce_write) then 654 call outwf(tdks%cg,dtset,psps,tdks%eigen0,fname,tdks%hdr,tdks%kg,dtset%kptns, & 655 & dtset%mband,tdks%mcg,dtset%mkmem,mpi_enreg,dtset%mpw,dtset%natom, & 656 & dtset%nband,dtset%nkpt,tdks%npwarr,dtset%nsppol,tdks%occ0,response, & 657 & dtfil%unwff2,tdks%wvl%wfs,tdks%wvl%descr, force_write=.TRUE.) 658 else 659 call outwf(tdks%cg,dtset,psps,tdks%eigen0,fname,tdks%hdr,tdks%kg,dtset%kptns, & 660 & dtset%mband,tdks%mcg,dtset%mkmem,mpi_enreg,dtset%mpw,dtset%natom, & 661 & dtset%nband,dtset%nkpt,tdks%npwarr,dtset%nsppol,tdks%occ0,response, & 662 & dtfil%unwff2,tdks%wvl%wfs,tdks%wvl%descr) 663 end if 664 665 end subroutine prt_wfk
m_rttddft_output/rttddft_output [ Functions ]
[ Top ] [ m_rttddft_output ] [ Functions ]
NAME
rttddft_output
FUNCTION
Main output subroutine
INPUTS
dtfil <type datafiles_type> = infos about file names, file unit numbers dtset <type(dataset_type)> = all input variables for this dataset istep <integer> = step number mpi_enreg <MPI_type> = MPI-parallelisation information psps <type(pseudopotential_type)> = variables related to pseudopotentials tdks <type(tdks_type)> = the tdks object to initialize
OUTPUT
SOURCE
80 subroutine rttddft_output(dtfil, dtset, istep, mpi_enreg, psps, tdks) 81 82 implicit none 83 84 !Arguments ------------------------------------ 85 !scalars 86 integer, intent(in) :: istep 87 type(datafiles_type), intent(inout) :: dtfil 88 type(dataset_type), intent(inout) :: dtset 89 type(MPI_type), intent(inout) :: mpi_enreg 90 type(pseudopotential_type), intent(inout) :: psps 91 type(tdks_type), intent(inout) :: tdks 92 !arrays 93 94 !Local variables------------------------------- 95 integer :: me 96 !scalars 97 character :: tmp 98 character(len=500) :: msg 99 integer :: stat 100 !arrays 101 102 ! ************************************************************************* 103 104 !** Special case of first step 105 if (istep == tdks%first_step) then 106 ! Open energy file and writes header if needed 107 if (open_file(tdks%fname_tdener, msg, newunit=tdks%tdener_unit, status='unknown', form='formatted') /= 0) then 108 write(msg,'(a,a)') 'Error while trying to open file ', tdks%fname_tdener 109 ABI_ERROR(msg) 110 end if 111 if (mpi_enreg%me == 0) then 112 if (dtset%td_restart>0) then 113 do 114 read(tdks%tdener_unit,*,iostat=stat) tmp 115 if (stat /= 0) exit 116 end do 117 backspace(tdks%tdener_unit) 118 else 119 write(msg,'(a)') "# RT-TDDFT -- Energy file. All quantities are in Hartree atomic units." 120 call wrtout(tdks%tdener_unit,msg) 121 write(msg,'(a)') "# step time E_total E_kinetic E_hartree E_xc E_ewald & 122 & E_corepsp E_localpsp E_nonlocalpsp E_paw E_entropy E_vdw" 123 call wrtout(tdks%tdener_unit,msg) 124 end if 125 end if 126 end if 127 128 !!FB: This is most probably not needed 129 !!Update header, with evolving variables 130 !call tdks%hdr%update(tdks%bantot,tdks%etot,tdks%energies%e_fermie,tdks%energies%e_fermih, & 131 ! & tdks%hdr%residm,tdks%rprimd,tdks%occ0,tdks%pawrhoij, & 132 ! & tdks%xred,dtset%amu_orig,comm_atom=mpi_enreg%comm_atom, & 133 ! & mpi_atmtab=mpi_enreg%my_atmtab) 134 135 !** Writes some info in main output file 136 write(msg,'(a,a,f14.6,a)') ch10, 'Total energy = ', tdks%etot,' Ha' 137 call wrtout(ab_out,msg) 138 if (do_write_log) call wrtout(std_out,msg) 139 140 write(msg,'(a,f14.6,a)') 'Integrated density (ie. total nb of electrons) = ', & 141 & sum(tdks%rhor(:,1))*tdks%ucvol/tdks%nfftf, ch10 142 call wrtout(ab_out,msg) 143 if (do_write_log) call wrtout(std_out,msg) 144 145 !** Writes in energy file 146 write(msg,'(i0,f10.5,11(f14.8,1X))') istep-1, (istep-1)*tdks%dt, tdks%etot, tdks%energies%e_kinetic, & 147 & tdks%energies%e_hartree, tdks%energies%e_xc, tdks%energies%e_ewald, & 148 & tdks%energies%e_corepsp, tdks%energies%e_localpsp, tdks%energies%e_nlpsp_vfock, & 149 & tdks%energies%e_paw, tdks%energies%e_entropy, tdks%energies%e_vdw_dftd 150 call wrtout(tdks%tdener_unit,msg) 151 152 !** Writes additional optional properties 153 !Computed at actual step 154 if (mod(istep,dtset%td_prtstr) == 0) then 155 call prt_den(dtfil,dtset,istep,mpi_enreg,psps,tdks) 156 end if 157 !Computed at previous step 158 if (mod(istep-1,dtset%td_prtstr) == 0) then 159 call prt_eig(dtfil,dtset,istep-1,mpi_enreg,tdks) 160 call prt_occ(dtfil,dtset,istep-1,mpi_enreg,tdks) 161 call prt_dos(dtfil,dtset,istep-1,mpi_enreg,psps,tdks) 162 end if 163 if (mod(istep,dtset%td_prtstr) == 0) then 164 if (dtset%prtwf > 0) then 165 call prt_wfk(dtfil,dtset,istep,mpi_enreg,psps,tdks) 166 call prt_restart(dtfil,istep,mpi_enreg,tdks) 167 end if 168 end if 169 170 !** Special case of last step 171 me = xmpi_comm_rank(mpi_enreg%comm_world) 172 if (istep == tdks%first_step+tdks%ntime-1) then 173 if (mod(istep,dtset%td_prtstr) /= 0 .or. dtset%prtwf <= 0) then 174 call prt_wfk(dtfil,dtset,istep,mpi_enreg,psps,tdks,force_write=.TRUE.) 175 call prt_restart(dtfil,istep,mpi_enreg,tdks) 176 end if 177 close(tdks%tdener_unit) 178 end if 179 180 end subroutine rttddft_output