Update WFN_MIX and RTP

src/cp_control_types.F

@@ -72,6 +72,7 @@ MODULE cp_control_types
       LOGICAL                              :: hfx_redistribute = .FALSE.
       INTEGER                              :: aspc_order = 0
       INTEGER                              :: sc_check_start = 0
+      LOGICAL                              :: apply_wfn_mix_init_restart = .FALSE.
       LOGICAL                              :: apply_delta_pulse = .FALSE.
       LOGICAL                              :: apply_delta_pulse_mag = .FALSE.
       LOGICAL                              :: periodic = .FALSE.

src/cp_control_utils.F

@@ -2262,6 +2262,8 @@ SUBROUTINE read_rtp_section(dft_control, rtp_section)
                                 i_val=dft_control%rtp_control%initial_wfn)
       CALL section_vals_val_get(rtp_section, "HFX_BALANCE_IN_CORE", &
                                 l_val=dft_control%rtp_control%hfx_redistribute)
+      CALL section_vals_val_get(rtp_section, "APPLY_WFN_MIX_INIT_RESTART", &
+                                l_val=dft_control%rtp_control%apply_wfn_mix_init_restart)
       CALL section_vals_val_get(rtp_section, "APPLY_DELTA_PULSE", &
                                 l_val=dft_control%rtp_control%apply_delta_pulse)
       CALL section_vals_val_get(rtp_section, "APPLY_DELTA_PULSE_MAG", &

src/emd/rt_propagation_utils.F

@@ -26,6 +26,7 @@ MODULE rt_propagation_utils
                                               cp_fm_to_fm,&
                                               cp_fm_type
    USE cp_log_handling,                 ONLY: cp_get_default_logger,&
+                                              cp_logger_get_default_io_unit,&
                                               cp_logger_get_default_unit_nr,&
                                               cp_logger_type
    USE cp_output_handling,              ONLY: cp_p_file,&
@@ -91,6 +92,7 @@ MODULE rt_propagation_utils
    USE qs_rho_types,                    ONLY: qs_rho_get,&
                                               qs_rho_set,&
                                               qs_rho_type
+   USE qs_scf_wfn_mix,                  ONLY: wfn_mix
    USE qs_subsys_types,                 ONLY: qs_subsys_get,&
                                               qs_subsys_type
    USE rt_propagation_types,            ONLY: get_rtp,&
@@ -261,8 +263,8 @@ SUBROUTINE get_restart_wfn(qs_env)
       TYPE(qs_environment_type), POINTER                 :: qs_env
 
       CHARACTER(LEN=default_path_length)                 :: file_name, project_name
-      INTEGER                                            :: i, id_nr, im, ispin, ncol, nspin, re, &
-                                                            unit_nr
+      INTEGER                                            :: i, id_nr, im, ispin, ncol, nspin, &
+                                                            output_unit, re, unit_nr
       REAL(KIND=dp)                                      :: alpha, cs_pos
       TYPE(atomic_kind_type), DIMENSION(:), POINTER      :: atomic_kind_set
       TYPE(cp_fm_type)                                   :: mos_occ
@@ -292,6 +294,8 @@ SUBROUTINE get_restart_wfn(qs_env)
                       rho=rho_struct, &
                       para_env=para_env)
       logger => cp_get_default_logger()
+      output_unit = cp_logger_get_default_io_unit(logger)
+
       IF (logger%para_env%is_source()) THEN
          unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.)
       ELSE
@@ -308,6 +312,10 @@ SUBROUTINE get_restart_wfn(qs_env)
                                        id_nr=id_nr, multiplicity=dft_control%multiplicity, dft_section=dft_section)
          CALL set_uniform_occupation_mo_array(mo_array, nspin)
 
+         IF (dft_control%rtp_control%apply_wfn_mix_init_restart) &
+            CALL wfn_mix(mo_array, particle_set, dft_section, qs_kind_set, para_env, output_unit, &
+                         for_rtp=.TRUE.)
+
          DO ispin = 1, nspin
             CALL calculate_density_matrix(mo_array(ispin), p_rmpv(ispin)%matrix)
          END DO

src/input_constants.F

@@ -926,6 +926,7 @@ MODULE input_constants
                                                ramp_env = 3, &
                                                custom_env = 4
 
+   ! type of reference for RTP%PRINT%PROJECTION_MO
    INTEGER, PARAMETER, PUBLIC               :: proj_mo_ref_scf = 1, &
                                                proj_mo_ref_xas_tdp = 2
 
@@ -1261,4 +1262,9 @@ MODULE input_constants
                                                scan_yz = 23, &
                                                scan_xyz = 123
 
+   ! DFT%PRINT%WFN_MIX orgine type
+   INTEGER, PARAMETER, PUBLIC               :: wfn_mix_orig_occ = 1, &
+                                               wfn_mix_orig_virtual = 2, &
+                                               wfn_mix_orig_external = 3
+
 END MODULE input_constants

src/input_cp2k_dft.F

@@ -101,10 +101,11 @@ MODULE input_cp2k_dft
         weight_type_unit, wfi_aspc_nr, wfi_frozen_method_nr, wfi_linear_p_method_nr, &
         wfi_linear_ps_method_nr, wfi_linear_wf_method_nr, wfi_ps_method_nr, &
         wfi_use_guess_method_nr, wfi_use_prev_p_method_nr, wfi_use_prev_rho_r_method_nr, &
-        wfi_use_prev_wf_method_nr, xas_1s_type, xas_2p_type, xas_2s_type, xas_3d_type, &
-        xas_3p_type, xas_3s_type, xas_4d_type, xas_4f_type, xas_4p_type, xas_4s_type, xas_dip_len, &
-        xas_dip_vel, xas_dscf, xas_none, xas_not_excited, xas_tdp_by_index, xas_tdp_by_kind, &
-        xas_tp_fh, xas_tp_flex, xas_tp_hh, xas_tp_xfh, xas_tp_xhh, xes_tp_val
+        wfi_use_prev_wf_method_nr, wfn_mix_orig_external, wfn_mix_orig_occ, wfn_mix_orig_virtual, &
+        xas_1s_type, xas_2p_type, xas_2s_type, xas_3d_type, xas_3p_type, xas_3s_type, xas_4d_type, &
+        xas_4f_type, xas_4p_type, xas_4s_type, xas_dip_len, xas_dip_vel, xas_dscf, xas_none, &
+        xas_not_excited, xas_tdp_by_index, xas_tdp_by_kind, xas_tp_fh, xas_tp_flex, xas_tp_hh, &
+        xas_tp_xfh, xas_tp_xhh, xes_tp_val
    USE input_cp2k_almo,                 ONLY: create_almo_scf_section
    USE input_cp2k_distribution,         ONLY: create_distribution_section
    USE input_cp2k_ec,                   ONLY: create_ec_section
@@ -2849,12 +2850,19 @@ SUBROUTINE create_wfn_mix_section(section)
                           description="A section that allows manipulation of the MO coeffs,"// &
                           " e.g. for changing a ground state into an excited state."// &
                           "Starting from a copy of the original MOs, changes can be made"// &
-                          "by adding linear combinations of HOMO/LUMO of the original MOs to the result MOs", &
+                          "by adding linear combinations of HOMO/LUMO of the original MOs to the result MOs."// &
+                          "This method is called after an SCF optimization or before an RTP run if "// &
+                          "INITIAL_WFN=RESTART_WFN. Note that if called after an SCF optimization, a restart file "// &
+                          "with the mixed MOs is saved. This is not the case for an RTP with "// &
+                          "INITIAL_WFN=RESTART_WFN.", &
                           n_keywords=1, n_subsections=0, repeats=.FALSE.)
 
       CALL keyword_create(keyword, __LOCATION__, name="OVERWRITE_MOS", &
-                          description="If the keyword is active molecular orbitals in memory will be replaced by the mixed wfn."// &
-                          " In combination with RTP or EMD no restart will be required to use the mixed wfn.", &
+                        description="If set to True, the active molecular orbitals in memory will be replaced by the mixed wfn "// &
+                        "at the end of the wfn mixing procedure. For instance, you can then use this new set of MOs to perform "// &
+                          "RTP or EMD directly. Not that in the case of an RTP run with INITIAL_WFN=RESTART_WFN, this keyword "// &
+                          "is not used: the mixed wfn is used to start the time-dependent run. "// &
+                          "Default value is False.", &
                           default_l_val=.FALSE., lone_keyword_l_val=.TRUE.)
       CALL section_add_keyword(section, keyword)
       CALL keyword_release(keyword)
@@ -2865,7 +2873,8 @@ SUBROUTINE create_wfn_mix_section(section)
                           n_keywords=1, n_subsections=0, repeats=.TRUE.)
 
       CALL keyword_create(keyword, __LOCATION__, name="RESULT_MO_INDEX", &
-                          description="Index of the MO (y) to be modified. Counting down in energy with HOMO=1", &
+                          description="Index of the MO (y) to be modified. Counting down in energy: "// &
+                          "set to 1 for the HOMO, to 3 for the HOMO-3.", &
                           usage="RESULT_MO_INDEX 1", type_of_var=integer_t, default_i_val=0)
       CALL section_add_keyword(subsection, keyword)
       CALL keyword_release(keyword)
@@ -2895,7 +2904,13 @@ SUBROUTINE create_wfn_mix_section(section)
 
       CALL keyword_create(keyword, __LOCATION__, name="ORIG_MO_INDEX", &
                           description="Index of the original MO (x). "// &
-                          "Counting down in energy with HOMO=1 or up from LUMO=1, depending on ORIG_IS_VIRTUAL.", &
+                          "If ORIG_TYPE is OCCUPIED, it counts down in energy: set to 1 to point to "// &
+                          "the HOMO and to 3 for the HOMO-3. "// &
+                          "If ORIG_TYPE is VIRTUAL, it counts up in energy: set to 1 to point to "// &
+                          "the LUMO and to 3 for the LUMO+3. "// &
+                          "If ORIG_TYPE is EXTERNAL, it counts down in energy for the external "// &
+                          "set of MOs: set to 1 to point to the HOMO and to 3 for the HOMO-3."// &
+                          "Do not set to zero or negative values.", &
                           usage="ORIG_MO_INDEX 1", type_of_var=integer_t, default_i_val=0)
       CALL section_add_keyword(subsection, keyword)
       CALL keyword_release(keyword)
@@ -2923,9 +2938,31 @@ SUBROUTINE create_wfn_mix_section(section)
       CALL section_add_keyword(subsection, keyword)
       CALL keyword_release(keyword)
 
-      CALL keyword_create(keyword, __LOCATION__, name="ORIG_IS_VIRTUAL", &
-                          description="The original MO (x) is a LUMO.", &
-                          usage="ORIG_IS_VIRTUAL", type_of_var=logical_t)
+      CALL keyword_create(keyword, __LOCATION__, name="ORIG_TYPE", &
+                          description="Type of the orgine MO. Note that using VIRTUAL along with an RTP run with "// &
+                          "INITIAL_WFN=RESTART_WFN will probably crash since no virtual state "// &
+                          "will be available.", &
+                          enum_c_vals=s2a("OCCUPIED", "VIRTUAL", 'EXTERNAL'), &
+                          usage="ORIG_TYPE OCCUPIED", &
+                          default_i_val=wfn_mix_orig_occ, &
+                          enum_desc=s2a("The original MO is in the MOs matrix itself.", &
+                                        "The original MO is not in the MOs matrix but in the 'unoccupied' one. "// &
+                                        "The SCF cycle shall provide it.", &
+                                        "The orginal MO is from an external .wfn file. Use the keyword "// &
+                                        "ORIG_EXT_FILE_NAME to define its name."), &
+                          enum_i_vals=(/wfn_mix_orig_occ, wfn_mix_orig_virtual, wfn_mix_orig_external/))
+      CALL section_add_keyword(subsection, keyword)
+      CALL keyword_release(keyword)
+
+      CALL keyword_create(keyword, __LOCATION__, name="ORIG_EXT_FILE_NAME", &
+                          description="Name of the wavefunction file to read the original MO from. "// &
+                          "For instance, a restart wfn file from SCF calculation or an excited state from XAS_TDP calculation. "// &
+                          "If no file is specified, the run will crash."// &
+                          "Currently, a RTP restart file (.rtpwfn) cannot be used as reference."// &
+                          "Currently, this file SHALL have the basis set, number of MO and the same number of spin as the one "// &
+                          "from the SCF cycle.", &
+                          usage="ORIG_EXT_FILE_NAME <FILENAME>", &
+                          default_lc_val="EMPTY")
       CALL section_add_keyword(subsection, keyword)
       CALL keyword_release(keyword)
 
@@ -8531,6 +8568,18 @@ SUBROUTINE create_rtp_section(section)
       CALL section_add_keyword(section, keyword)
       CALL keyword_release(keyword)
 
+      CALL keyword_create(keyword, __LOCATION__, name="APPLY_WFN_MIX_INIT_RESTART", &
+                          description="If set to True and in the case of INITIAL_WFN=RESTART_WFN, call the "// &
+                          "DFT%PRINT%WFN_MIX section to mix the read initial wfn. The starting wave-function of the "// &
+                          "RTP will be the mixed one. Setting this to True without a defined WFN_MIX section will "// &
+                          "not do anything as defining a WFN_MIX section without this keyword for RTP run with "// &
+                          "INITIAL_WFN=RESTART_WFN. Note that if INITIAL_WFN=SCF_WFN, this keyword is not needed to "// &
+                          "apply the mixing defined in the WFN_MIX section. Default is False.", &
+                          usage="APPLY_WFN_MIX_INIT_RESTART", &
+                          default_l_val=.FALSE., lone_keyword_l_val=.TRUE.)
+      CALL section_add_keyword(section, keyword)
+      CALL keyword_release(keyword)
+
       CALL keyword_create(keyword, __LOCATION__, name="APPLY_DELTA_PULSE", &
                           description="Applies a delta kick to the initial wfn (only RTP for now - the EMD "// &
                           " case is not yet implemented). Only work for INITIAL_WFN=SCF_WFN", &

src/input_cp2k_projection_rtp.F

@@ -47,10 +47,12 @@ SUBROUTINE create_projection_rtp_section(section)
 
       CPASSERT(.NOT. ASSOCIATED(section))
       CALL section_create(section, __LOCATION__, name="PROJECTION_MO", &
-                          description="Projects the time dependent MO "// &
-                          "coefficients on reference ones. You can define "// &
+                          description="Projects the Time Dependent (TD) MO "// &
+                          "coefficients to reference ones. You can define "// &
                           "several sections like this to project the TD-MOs "// &
-                          "on different reference MOs.", &
+                          "on different reference MOs. Note that each section "// &
+                          "projects from one spin of the TD MOs (TD_MO_INDEX) to "// &
+                          "one spin of the reference MOs (REF_MO_INDEX).", &
                           n_keywords=7, n_subsections=1, &
                           repeats=.TRUE.)
 
@@ -67,10 +69,12 @@ SUBROUTINE create_projection_rtp_section(section)
       CALL keyword_release(keyword)
 
       CALL keyword_create(keyword, __LOCATION__, name="REF_MO_FILE_NAME", &
-                          description="Name of the wavefunction file to read the reference MO form. "// &
-                          "For instance a restart wfn file from SCF calculation or an excited state from XAS_TDP calculation. "// &
+                          description="Name of the wavefunction file to read the reference MO from. "// &
+                          "For instance, a restart wfn file from SCF calculation or an excited state from XAS_TDP calculation. "// &
                           "If no file is specified, the default is to use DFT%WFN_RESTART_FILE_NAME."// &
-                          "Currently, a RTP restart file (.rtpwfn) cannot be used as reference.", &
+                          "Currently, a RTP restart file (.rtpwfn) cannot be used as reference."// &
+                          "Currently, this file SHALL have the same number of spin as the propagated one "// &
+                          "(eventhough you use only the first spin from this reference).", &
                           usage="REF_MO_FILE_NAME <FILENAME>", &
                           default_lc_val="DEFAULT")
       CALL section_add_keyword(section, keyword)
@@ -91,7 +95,7 @@ SUBROUTINE create_projection_rtp_section(section)
                           description="Spin of the reference MOs to consider. "// &
                           "1 for ALPHA and 2 for BETA spin respectively. "// &
                           "If the reference MO is spin independent this key is not used.", &
-                          usage="MO_TD_SPIN 1", &
+                          usage="REF_MO_SPIN 1", &
                           default_i_val=1, &
                           n_var=1, type_of_var=integer_t, repeats=.FALSE.)
 
@@ -119,7 +123,7 @@ SUBROUTINE create_projection_rtp_section(section)
       CALL keyword_create(keyword, __LOCATION__, name="TD_MO_SPIN", &
                           description="Spin of the TD MOs to consider. 1 for ALPHA spin, 2 for BETA spin. "// &
                           "If the TD calculation is spin independent this key is not used.", &
-                          usage="MO_TD_SPIN 1 2", &
+                          usage="MO_TD_SPIN 1", &
                           default_i_val=1, &
                           n_var=1, type_of_var=integer_t)
       CALL section_add_keyword(section, keyword)

src/qs_scf_post_gpw.F

@@ -272,6 +272,7 @@ SUBROUTINE scf_post_calculation_gpw(qs_env, wf_type, do_mp2)
       TYPE(dft_control_type), POINTER                    :: dft_control
       TYPE(mo_set_type), DIMENSION(:), POINTER           :: mos
       TYPE(molecule_type), POINTER                       :: molecule_set(:)
+      TYPE(mp_para_env_type), POINTER                    :: para_env
       TYPE(particle_list_type), POINTER                  :: particles
       TYPE(particle_type), DIMENSION(:), POINTER         :: particle_set
       TYPE(pw_env_type), POINTER                         :: pw_env
@@ -749,9 +750,11 @@ SUBROUTINE scf_post_calculation_gpw(qs_env, wf_type, do_mp2)
       IF (do_kpoints) THEN
          ! nothing at the moment, not implemented
       ELSE
-         CALL get_qs_env(qs_env, matrix_s=matrix_s)
-         CALL wfn_mix(mos, particle_set, dft_section, qs_kind_set, &
-                      lumo_ptr, scf_env, matrix_s, output_unit, marked_states)
+         CALL get_qs_env(qs_env, matrix_s=matrix_s, para_env=para_env)
+         CALL wfn_mix(mos, particle_set, dft_section, qs_kind_set, para_env, &
+                      output_unit, unoccupied_orbs=lumo_ptr, scf_env=scf_env, &
+                      matrix_s=matrix_s, marked_states=marked_states)
+
          IF (p_loc_lumo) CALL cp_fm_release(lumo_localized)
       END IF
       IF (ASSOCIATED(marked_states)) THEN