Skip to content

Commit

Permalink
XAS_TDP| Allow the dumping of an excited state AMEW into a restart file
Browse files Browse the repository at this point in the history
  • Loading branch information
@abussy
abussy committed Jul 29, 2022
1 parent 4a74572 commit 840caac
Show file tree
Hide file tree
Showing 2 changed files with 107 additions and 19 deletions.
14 changes: 14 additions & 0 deletions src/input_cp2k_dft.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -7973,6 +7973,20 @@ SUBROUTINE create_xas_tdp_section(section)
CALL section_add_subsection(subsection, print_key)
CALL section_release(print_key)

CALL cp_print_key_section_create(print_key, __LOCATION__, name="RESTART_WFN", &
description="Controles the dumping of a MO restart file for a given "// &
"excited state index. Only for K-edge RKS calculations. "// &
"Can be repeated to get multiple *.wfn files at once.", &
print_level=debug_print_level, filename="", common_iter_levels=1)
CALL keyword_create(keyword, __LOCATION__, name="EXCITED_STATE_INDEX", variants=(/"INDEX"/), &
description="The index of the excited state that should be dumped", &
usage="INDEX {int}", default_i_val=1, repeats=.TRUE.)
CALL section_add_keyword(print_key, keyword)
CALL keyword_release(keyword)

CALL section_add_subsection(subsection, print_key)
CALL section_release(print_key)

CALL create_pdos_section(print_key)
CALL section_add_subsection(subsection, print_key)
CALL section_release(print_key)
Expand Down
112 changes: 93 additions & 19 deletions src/xas_tdp_methods.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -104,9 +104,11 @@ MODULE xas_tdp_methods
USE qs_loc_utils, ONLY: qs_loc_control_init,&
qs_loc_env_init,&
set_loc_centers
USE qs_mo_io, ONLY: write_mo_set_low
USE qs_mo_methods, ONLY: calculate_subspace_eigenvalues
USE qs_mo_types, ONLY: allocate_mo_set,&
deallocate_mo_set,&
duplicate_mo_set,&
get_mo_set,&
init_mo_set,&
mo_set_p_type,&
Expand Down Expand Up @@ -453,7 +455,7 @@ SUBROUTINE xas_tdp_core(xas_tdp_section, qs_env)
xas_tdp_control, xas_tdp_env)
CALL xas_tdp_post(tddfpt_spin_cons, current_state, xas_tdp_env, &
xas_tdp_section, qs_env)
CALL write_donor_state_restart(tddfpt_spin_cons, current_state, xas_tdp_section)
CALL write_donor_state_restart(tddfpt_spin_cons, current_state, xas_tdp_section, qs_env)
END IF
IF (xas_tdp_control%do_spin_flip) THEN
Expand All @@ -462,7 +464,7 @@ SUBROUTINE xas_tdp_core(xas_tdp_section, qs_env)
!no dipole in spin-flip (spin-forbidden)
CALL xas_tdp_post(tddfpt_spin_flip, current_state, xas_tdp_env, &
xas_tdp_section, qs_env)
CALL write_donor_state_restart(tddfpt_spin_flip, current_state, xas_tdp_section)
CALL write_donor_state_restart(tddfpt_spin_flip, current_state, xas_tdp_section, qs_env)
END IF
IF (xas_tdp_control%do_singlet) THEN
Expand All @@ -473,7 +475,7 @@ SUBROUTINE xas_tdp_core(xas_tdp_section, qs_env)
xas_tdp_control, xas_tdp_env)
CALL xas_tdp_post(tddfpt_singlet, current_state, xas_tdp_env, &
xas_tdp_section, qs_env)
CALL write_donor_state_restart(tddfpt_singlet, current_state, xas_tdp_section)
CALL write_donor_state_restart(tddfpt_singlet, current_state, xas_tdp_section, qs_env)
END IF
IF (xas_tdp_control%do_triplet) THEN
Expand All @@ -482,7 +484,7 @@ SUBROUTINE xas_tdp_core(xas_tdp_section, qs_env)
!no dipole for triplets by construction
CALL xas_tdp_post(tddfpt_triplet, current_state, xas_tdp_env, &
xas_tdp_section, qs_env)
CALL write_donor_state_restart(tddfpt_triplet, current_state, xas_tdp_section)
CALL write_donor_state_restart(tddfpt_triplet, current_state, xas_tdp_section, qs_env)
END IF
! Include the SOC if required, only for 2p donor stataes
Expand Down Expand Up @@ -2444,12 +2446,11 @@ SUBROUTINE xas_tdp_post(ex_type, donor_state, xas_tdp_env, xas_tdp_section, qs_e
CHARACTER(len=*), PARAMETER :: routineN = 'xas_tdp_post'
CHARACTER(len=default_string_length) :: domo, domon, excite, pos, xas_mittle
INTEGER :: handle, ic, ido_mo, imo, irep, ispin, &
n_dependent, n_rep, nao, ncubes, &
ndo_mo, ndo_so, nlumo, nmo, nspins, &
output_unit
INTEGER :: ex_state_idx, handle, ic, ido_mo, imo, irep, ispin, n_dependent, n_rep, nao, &
ncubes, ndo_mo, ndo_so, nlumo, nmo, nspins, output_unit
INTEGER, DIMENSION(:), POINTER :: bounds, list, state_list
LOGICAL :: append_cube, do_cubes, do_pdos
LOGICAL :: append_cube, do_cubes, do_pdos, &
do_wfn_restart
REAL(dp), DIMENSION(:), POINTER :: lr_evals
REAL(dp), DIMENSION(:, :), POINTER :: centers
TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
Expand All @@ -2459,33 +2460,37 @@ SUBROUTINE xas_tdp_post(ex_type, donor_state, xas_tdp_env, xas_tdp_section, qs_e
TYPE(cp_logger_type), POINTER :: logger
TYPE(cp_para_env_type), POINTER :: para_env
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
TYPE(mo_set_p_type), DIMENSION(:), POINTER :: mos, restart_mos
TYPE(mo_set_type), POINTER :: mo_set
TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
TYPE(section_vals_type), POINTER :: print_key
NULLIFY (atomic_kind_set, particle_set, qs_kind_set, mo_set, lr_evals, lr_coeffs, mo_coeff)
NULLIFY (mo_struct, para_env, blacs_env, fm_struct, matrix_s, work_fm, print_key, logger)
NULLIFY (bounds, state_list, list)
NULLIFY (bounds, state_list, list, restart_mos, mos)
!Tests on what to do
logger => cp_get_default_logger()
do_pdos = .FALSE.; do_cubes = .FALSE.
do_pdos = .FALSE.; do_cubes = .FALSE.; do_wfn_restart = .FALSE.
IF (BTEST(cp_print_key_should_output(logger%iter_info, xas_tdp_section, &
"PRINT%PDOS"), cp_p_file)) do_pdos = .TRUE.
IF (BTEST(cp_print_key_should_output(logger%iter_info, xas_tdp_section, &
"PRINT%CUBES"), cp_p_file)) do_cubes = .TRUE.
IF (.NOT. (do_pdos .OR. do_cubes)) RETURN
IF (BTEST(cp_print_key_should_output(logger%iter_info, xas_tdp_section, &
"PRINT%RESTART_WFN"), cp_p_file)) do_wfn_restart = .TRUE.
IF (.NOT. (do_pdos .OR. do_cubes .OR. do_wfn_restart)) RETURN
CALL timeset(routineN, handle)
!Initialization
CALL get_qs_env(qs_env, atomic_kind_set=atomic_kind_set, particle_set=particle_set, &
qs_kind_set=qs_kind_set, para_env=para_env, blacs_env=blacs_env, &
matrix_s=matrix_s)
matrix_s=matrix_s, mos=mos)
SELECT CASE (ex_type)
CASE (tddfpt_spin_cons)
Expand Down Expand Up @@ -2528,6 +2533,47 @@ SUBROUTINE xas_tdp_post(ex_type, donor_state, xas_tdp_env, xas_tdp_section, qs_e
nrow_global=nao, ncol_global=nmo)
CALL cp_fm_create(mo_coeff, mo_struct)
!Dump the TDDFT excited state AMEW wavefunction into a file for restart in RTP
IF (do_wfn_restart) THEN
IF (.NOT. (nspins == 1 .AND. donor_state%state_type == xas_1s_type)) THEN
CPABORT("RESTART.wfn file only available for RKS K-edge XAS spectroscopy")
END IF
CALL section_vals_val_get(xas_tdp_section, "PRINT%RESTART_WFN%EXCITED_STATE_INDEX", n_rep_val=n_rep)
DO irep = 1, n_rep
CALL section_vals_val_get(xas_tdp_section, "PRINT%RESTART_WFN%EXCITED_STATE_INDEX", &
i_rep_val=irep, i_val=ex_state_idx)
CPASSERT(ex_state_idx <= SIZE(lr_evals))
ALLOCATE (restart_mos(2))
DO ispin = 1, 2
CALL duplicate_mo_set(restart_mos(ispin)%mo_set, mos(1)%mo_set)
END DO
CALL cp_fm_to_fm_submat(msource=lr_coeffs, mtarget=restart_mos(1)%mo_set%mo_coeff, nrow=nao, &
ncol=1, s_firstrow=1, s_firstcol=ex_state_idx, t_firstrow=1, &
t_firstcol=donor_state%mo_indices(1, 1))
xas_mittle = 'xasat'//TRIM(ADJUSTL(cp_to_string(donor_state%at_index)))//'_'//TRIM(domo)// &
'_'//TRIM(excite)//'_idx'//TRIM(ADJUSTL(cp_to_string(ex_state_idx)))
output_unit = cp_print_key_unit_nr(logger, xas_tdp_section, "PRINT%RESTART_WFN", &
extension=".wfn", file_status="REPLACE", &
file_action="WRITE", file_form="UNFORMATTED", &
middle_name=xas_mittle)
CALL write_mo_set_low(restart_mos, particle_set=particle_set, &
qs_kind_set=qs_kind_set, ires=output_unit)
CALL cp_print_key_finished_output(output_unit, logger, xas_tdp_section, "PRINT%RESTART_WFN")
DO ispin = 1, 2
CALL deallocate_mo_set(restart_mos(ispin)%mo_set)
END DO
DEALLOCATE (restart_mos)
END DO
END IF
!PDOS related stuff
IF (do_pdos) THEN
Expand Down Expand Up @@ -3149,27 +3195,31 @@ END SUBROUTINE print_xas_tdp_to_file
!> \param ex_type singlet, triplet, etc.
!> \param donor_state ...
!> \param xas_tdp_section ...
!> \param qs_env ...
! **************************************************************************************************
SUBROUTINE write_donor_state_restart(ex_type, donor_state, xas_tdp_section)
SUBROUTINE write_donor_state_restart(ex_type, donor_state, xas_tdp_section, qs_env)
INTEGER, INTENT(IN) :: ex_type
TYPE(donor_state_type), POINTER :: donor_state
TYPE(section_vals_type), POINTER :: xas_tdp_section
TYPE(qs_environment_type), POINTER :: qs_env
CHARACTER(len=*), PARAMETER :: routineN = 'write_donor_state_restart'
CHARACTER(len=default_path_length) :: filename
CHARACTER(len=default_string_length) :: domo, excite, my_middle
INTEGER :: ex_atom, handle, nao, ndo_mo, nex, &
nspins, output_unit, rst_unit, &
INTEGER :: ex_atom, handle, ispin, nao, ndo_mo, &
nex, nspins, output_unit, rst_unit, &
state_type
INTEGER, DIMENSION(:, :), POINTER :: mo_indices
LOGICAL :: do_print
REAL(dp), DIMENSION(:), POINTER :: lr_evals
TYPE(cp_fm_type), POINTER :: lr_coeffs
TYPE(cp_logger_type), POINTER :: logger
TYPE(mo_set_p_type), DIMENSION(:), POINTER :: mos
TYPE(section_vals_type), POINTER :: print_key
NULLIFY (logger, lr_coeffs, lr_evals, print_key)
NULLIFY (logger, lr_coeffs, lr_evals, print_key, mos)
!Initialization
logger => cp_get_default_logger()
Expand Down Expand Up @@ -3221,6 +3271,7 @@ SUBROUTINE write_donor_state_restart(ex_type, donor_state, xas_tdp_section)
CALL cp_fm_get_info(lr_coeffs, nrow_global=nao)
state_type = donor_state%state_type
ex_atom = donor_state%at_index
mo_indices => donor_state%mo_indices
!Opening restart file
rst_unit = -1
Expand All @@ -3242,10 +3293,17 @@ SUBROUTINE write_donor_state_restart(ex_type, donor_state, xas_tdp_section)
IF (rst_unit > 0) THEN
WRITE (rst_unit) ex_atom, state_type, ndo_mo, ex_type
WRITE (rst_unit) nao, nex, nspins
WRITE (rst_unit) mo_indices(:, :)
WRITE (rst_unit) lr_evals(:)
END IF
CALL cp_fm_write_unformatted(lr_coeffs, rst_unit)
!The MOs as well (because the may have been localized)
CALL get_qs_env(qs_env, mos=mos)
DO ispin = 1, nspins
CALL cp_fm_write_unformatted(mos(ispin)%mo_set%mo_coeff, rst_unit)
END DO
!closing
CALL cp_print_key_finished_output(rst_unit, logger, xas_tdp_section, "PRINT%RESTART")
Expand All @@ -3269,17 +3327,19 @@ SUBROUTINE read_donor_state_restart(donor_state, ex_type, filename, qs_env)
CHARACTER(len=*), PARAMETER :: routineN = 'read_donor_state_restart'
INTEGER :: group, handle, nao, nex, nspins, &
INTEGER :: group, handle, ispin, nao, nex, nspins, &
output_unit, read_params(7), rst_unit, &
source
INTEGER, DIMENSION(:, :), POINTER :: mo_indices
LOGICAL :: file_exists
REAL(dp), DIMENSION(:), POINTER :: lr_evals
TYPE(cp_blacs_env_type), POINTER :: blacs_env
TYPE(cp_fm_struct_type), POINTER :: fm_struct
TYPE(cp_fm_type), POINTER :: lr_coeffs
TYPE(cp_para_env_type), POINTER :: para_env
TYPE(mo_set_p_type), DIMENSION(:), POINTER :: mos
NULLIFY (lr_evals, lr_coeffs, para_env, fm_struct, blacs_env)
NULLIFY (lr_evals, lr_coeffs, para_env, fm_struct, blacs_env, mos)
CALL timeset(routineN, handle)
Expand Down Expand Up @@ -3321,6 +3381,13 @@ SUBROUTINE read_donor_state_restart(donor_state, ex_type, filename, qs_env)
nex = read_params(6)
nspins = read_params(7)
ALLOCATE (mo_indices(donor_state%ndo_mo, nspins))
IF (rst_unit > 0) THEN
READ (rst_unit) mo_indices(1:donor_state%ndo_mo, 1:nspins)
END IF
CALL mp_bcast(mo_indices, source, group)
donor_state%mo_indices => mo_indices
!read evals
ALLOCATE (lr_evals(nex))
IF (rst_unit > 0) READ (rst_unit) lr_evals(1:nex)
Expand All @@ -3333,6 +3400,12 @@ SUBROUTINE read_donor_state_restart(donor_state, ex_type, filename, qs_env)
CALL cp_fm_read_unformatted(lr_coeffs, rst_unit)
CALL cp_fm_struct_release(fm_struct)
!read MO coeffs and replace in qs_env
CALL get_qs_env(qs_env, mos=mos)
DO ispin = 1, nspins
CALL cp_fm_read_unformatted(mos(ispin)%mo_set%mo_coeff, rst_unit)
END DO
!closing file
IF (para_env%ionode) THEN
CALL close_file(unit_number=rst_unit)
Expand Down Expand Up @@ -3388,6 +3461,7 @@ SUBROUTINE restart_calculation(rst_filename, xas_tdp_section, qs_env)
!clean-up
CALL xas_tdp_env_release(xas_tdp_env)
CALL free_ds_memory(donor_state)
DEALLOCATE (donor_state%mo_indices)
DEALLOCATE (donor_state)
END SUBROUTINE restart_calculation
Expand Down

0 comments on commit 840caac

Please sign in to comment.