create TDDFPT%LINRES section (#3101)

Co-authored-by: BelizSertcan <beliz.sertcan@uzh.ch>
cp2k · Nov 10, 2023 · 14576fc · 14576fc
1 parent 2f34a03
commit 14576fc
Show file tree

Hide file tree

Showing 14 changed files with 155 additions and 56 deletions.
diff --git a/src/excited_states.F b/src/excited_states.F
@@ -20,6 +20,8 @@ MODULE excited_states
                                               cp_logger_type
    USE ex_property_calculation,         ONLY: ex_properties
    USE exstates_types,                  ONLY: excited_energy_type
+   USE input_section_types,             ONLY: section_vals_get_subs_vals,&
+                                              section_vals_type
    USE qs_energy_types,                 ONLY: qs_energy_type
    USE qs_environment_types,            ONLY: get_qs_env,&
                                               qs_environment_type,&
@@ -72,6 +74,7 @@ SUBROUTINE excited_state_energy(qs_env, calculate_forces)
       TYPE(qs_energy_type), POINTER                      :: energy
       TYPE(qs_force_type), DIMENSION(:), POINTER         :: ks_force, lr_force
       TYPE(qs_p_env_type)                                :: p_env
+      TYPE(section_vals_type), POINTER                   :: tdlr_section
 
       CALL timeset(routineN, handle)
 
@@ -102,7 +105,8 @@ SUBROUTINE excited_state_energy(qs_env, calculate_forces)
             CALL zero_qs_force(lr_force)
             CALL set_qs_env(qs_env, force=lr_force)
             !
-            CALL response_equation(qs_env, p_env, ex_env%cpmos, unit_nr)
+            tdlr_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%TDDFPT%LINRES")
+            CALL response_equation(qs_env, p_env, ex_env%cpmos, unit_nr, tdlr_section)
             !
             CALL get_qs_env(qs_env, dft_control=dft_control)
             IF (dft_control%qs_control%semi_empirical) THEN

diff --git a/src/input_cp2k_properties_dft.F b/src/input_cp2k_properties_dft.F
@@ -93,7 +93,7 @@ SUBROUTINE create_properties_section(section)
 
       NULLIFY (subsection, keyword)
 
-      CALL create_linres_section(subsection)
+      CALL create_linres_section(subsection, create_subsections=.TRUE.)
       CALL section_add_subsection(section, subsection)
       CALL section_release(subsection)
 
@@ -144,15 +144,35 @@ END SUBROUTINE create_properties_section
 !>      a linear response calculation
 !>      Available properties : none
 !> \param section the section to create
+!> \param create_subsections indicates whether or not subsections should be created
+!> \param default_set_tdlr default parameters to be used if called from TDDFPT
 !> \author MI
 ! **************************************************************************************************
-   SUBROUTINE create_linres_section(section)
+   SUBROUTINE create_linres_section(section, create_subsections, default_set_tdlr)
       TYPE(section_type), POINTER                        :: section
+      LOGICAL, INTENT(in)                                :: create_subsections
+      LOGICAL, INTENT(IN), OPTIONAL                      :: default_set_tdlr
 
+      INTEGER                                            :: def_max_iter, def_precond
+      REAL(KIND=DP)                                      :: def_egap, def_eps, def_eps_filter
       TYPE(keyword_type), POINTER                        :: keyword
       TYPE(section_type), POINTER                        :: print_key, subsection
 
-      NULLIFY (keyword, subsection, print_key)
+      NULLIFY (keyword, print_key)
+
+      IF (PRESENT(default_set_tdlr)) THEN
+         def_egap = 0.02_dp
+         def_eps = 1.0e-10_dp
+         def_eps_filter = 1.0e-15_dp
+         def_max_iter = 100
+         def_precond = ot_precond_full_single_inverse
+      ELSE
+         def_egap = 0.2_dp
+         def_eps = 1.e-6_dp
+         def_eps_filter = 0.0_dp
+         def_max_iter = 50
+         def_precond = ot_precond_none
+      END IF
 
       CPASSERT(.NOT. ASSOCIATED(section))
       CALL section_create(section, __LOCATION__, name="linres", &
@@ -163,19 +183,19 @@ SUBROUTINE create_linres_section(section)
 
       CALL keyword_create(keyword, __LOCATION__, name="EPS", &
                           description="target accuracy for the convergence of the conjugate gradient.", &
-                          usage="EPS 1.e-6", default_r_val=1.e-6_dp)
+                          usage="EPS 1.e-6", default_r_val=def_eps)
       CALL section_add_keyword(section, keyword)
       CALL keyword_release(keyword)
 
       CALL keyword_create(keyword, __LOCATION__, name="EPS_FILTER", &
                           description="Filter threshold for response density matrix.", &
-                          usage="EPS 1.e-8", default_r_val=0.0_dp)
+                          usage="EPS 1.e-8", default_r_val=def_eps_filter)
       CALL section_add_keyword(section, keyword)
       CALL keyword_release(keyword)
 
       CALL keyword_create(keyword, __LOCATION__, name="MAX_ITER", &
                           description="Maximum number of conjugate gradient iteration to be performed for one optimization.", &
-                          usage="MAX_ITER 200", default_i_val=50)
+                          usage="MAX_ITER 200", default_i_val=def_max_iter)
       CALL section_add_keyword(section, keyword)
       CALL keyword_release(keyword)
 
@@ -191,7 +211,7 @@ SUBROUTINE create_linres_section(section)
          "They differ in effectiveness, cost of construction, cost of application. "// &
          "Properly preconditioned minimization can be orders of magnitude faster than doing nothing.", &
          usage="PRECONDITIONER FULL_ALL", &
-         default_i_val=ot_precond_none, &
+         default_i_val=def_precond, &
          enum_c_vals=s2a("FULL_ALL", "FULL_SINGLE_INVERSE", "FULL_SINGLE", "FULL_KINETIC", "FULL_S_INVERSE", &
                          "NONE"), &
          enum_desc=s2a("Most effective state selective preconditioner based on diagonalization, "// &
@@ -214,7 +234,7 @@ SUBROUTINE create_linres_section(section)
       CALL keyword_create(keyword, __LOCATION__, name="ENERGY_GAP", &
                           description="Energy gap estimate [a.u.] for preconditioning", &
                           usage="ENERGY_GAP 0.1", &
-                          default_r_val=0.2_dp)
+                          default_r_val=def_egap)
       CALL section_add_keyword(section, keyword)
       CALL keyword_release(keyword)
 
@@ -234,62 +254,67 @@ SUBROUTINE create_linres_section(section)
       CALL section_add_keyword(section, keyword)
       CALL keyword_release(keyword)
 
-      CALL create_localize_section(subsection)
-      CALL section_add_subsection(section, subsection)
-      CALL section_release(subsection)
+      IF (create_subsections) THEN
+         NULLIFY (subsection)
 
-      CALL create_current_section(subsection)
-      CALL section_add_subsection(section, subsection)
-      CALL section_release(subsection)
+         CALL create_localize_section(subsection)
+         CALL section_add_subsection(section, subsection)
+         CALL section_release(subsection)
 
-      CALL create_nmr_section(subsection)
-      CALL section_add_subsection(section, subsection)
-      CALL section_release(subsection)
+         CALL create_current_section(subsection)
+         CALL section_add_subsection(section, subsection)
+         CALL section_release(subsection)
 
-      CALL create_spin_spin_section(subsection)
-      CALL section_add_subsection(section, subsection)
-      CALL section_release(subsection)
+         CALL create_nmr_section(subsection)
+         CALL section_add_subsection(section, subsection)
+         CALL section_release(subsection)
 
-      CALL create_epr_section(subsection)
-      CALL section_add_subsection(section, subsection)
-      CALL section_release(subsection)
+         CALL create_spin_spin_section(subsection)
+         CALL section_add_subsection(section, subsection)
+         CALL section_release(subsection)
 
-      CALL create_polarizability_section(subsection)
-      CALL section_add_subsection(section, subsection)
-      CALL section_release(subsection)
+         CALL create_epr_section(subsection)
+         CALL section_add_subsection(section, subsection)
+         CALL section_release(subsection)
 
-      CALL create_dcdr_section(subsection)
-      CALL section_add_subsection(section, subsection)
-      CALL section_release(subsection)
+         CALL create_polarizability_section(subsection)
+         CALL section_add_subsection(section, subsection)
+         CALL section_release(subsection)
 
-      CALL create_vcd_section(subsection)
-      CALL section_add_subsection(section, subsection)
-      CALL section_release(subsection)
+         CALL create_dcdr_section(subsection)
+         CALL section_add_subsection(section, subsection)
+         CALL section_release(subsection)
 
-      CALL section_create(subsection, __LOCATION__, name="PRINT", &
-                          description="printing of information during the linear response calculation", &
-                          repeats=.FALSE.)
+         CALL create_vcd_section(subsection)
+         CALL section_add_subsection(section, subsection)
+         CALL section_release(subsection)
 
-      CALL cp_print_key_section_create( &
-         print_key, __LOCATION__, "program_run_info", &
-         description="Controls the printing of basic iteration information during the LINRES calculation", &
-         print_level=low_print_level, add_last=add_last_numeric, filename="__STD_OUT__")
-      CALL section_add_subsection(subsection, print_key)
-      CALL section_release(print_key)
+         CALL section_create(subsection, __LOCATION__, name="PRINT", &
+                             description="printing of information during the linear response calculation", &
+                             repeats=.FALSE.)
 
-      CALL cp_print_key_section_create(print_key, __LOCATION__, "RESTART", &
-                                       description="Controls the dumping of restart file of the response wavefunction. "// &
-                                       "For each set of response functions, i.e. for each perturbation, "// &
-                                       "one different restart file is dumped. These restart files should be "// &
-                                       "employed only to restart the same type of LINRES calculation, "// &
-                                       "i.e. with the same perturbation.", &
-                                       print_level=low_print_level, common_iter_levels=3, each_iter_names=s2a("ITER"), &
-                                       add_last=add_last_numeric, each_iter_values=(/3/), filename="")
-      CALL section_add_subsection(subsection, print_key)
-      CALL section_release(print_key)
+         CALL cp_print_key_section_create( &
+            print_key, __LOCATION__, "program_run_info", &
+            description="Controls the printing of basic iteration information during the LINRES calculation", &
+            print_level=low_print_level, add_last=add_last_numeric, filename="__STD_OUT__")
+         CALL section_add_subsection(subsection, print_key)
+         CALL section_release(print_key)
 
-      CALL section_add_subsection(section, subsection)
-      CALL section_release(subsection)
+         CALL cp_print_key_section_create(print_key, __LOCATION__, "RESTART", &
+                                          description="Controls the dumping of restart file of the response wavefunction. "// &
+                                          "For each set of response functions, i.e. for each perturbation, "// &
+                                          "one different restart file is dumped. These restart files should be "// &
+                                          "employed only to restart the same type of LINRES calculation, "// &
+                                          "i.e. with the same perturbation.", &
+                                          print_level=low_print_level, common_iter_levels=3, each_iter_names=s2a("ITER"), &
+                                          add_last=add_last_numeric, each_iter_values=(/3/), filename="")
+         CALL section_add_subsection(subsection, print_key)
+         CALL section_release(print_key)
+
+         CALL section_add_subsection(section, subsection)
+         CALL section_release(subsection)
+
+      END IF
 
    END SUBROUTINE create_linres_section
 
@@ -1656,6 +1681,11 @@ SUBROUTINE create_tddfpt2_section(section)
       CALL section_add_subsection(section, subsection)
       CALL section_release(subsection)
 
+      ! LINRES section
+      CALL create_linres_section(subsection, create_subsections=.FALSE., default_set_tdlr=.TRUE.)
+      CALL section_add_subsection(section, subsection)
+      CALL section_release(subsection)
+
       ! PRINT subsection
       CALL section_create(subsection, __LOCATION__, name="PRINT", &
                           description="Printing of information during the TDDFT run.", repeats=.FALSE.)

diff --git a/src/qs_tddfpt2_fprint.F b/src/qs_tddfpt2_fprint.F
@@ -124,7 +124,7 @@ SUBROUTINE tddfpt_print_forces(qs_env, evects, evals, ostrength, print_section,
       TYPE(mp_para_env_type), POINTER                    :: para_env
       TYPE(qs_force_type), DIMENSION(:), POINTER         :: gs_force, ks_force, td_force
       TYPE(qs_p_env_type)                                :: p_env
-      TYPE(section_vals_type), POINTER                   :: force_section
+      TYPE(section_vals_type), POINTER                   :: force_section, tdlr_section
 
       CALL timeset(routineN, handle)
 
@@ -149,6 +149,7 @@ SUBROUTINE tddfpt_print_forces(qs_env, evects, evals, ostrength, print_section,
          END IF
          force_section => section_vals_get_subs_vals(print_section, "FORCES")
          CALL section_vals_val_get(force_section, "THRESHOLD", r_val=threshold)
+         tdlr_section => section_vals_get_subs_vals(qs_env%input, "PROPERTIES%TDDFPT%LINRES")
          ALLOCATE (state_list(nstates))
          CALL build_state_list(force_section, state_list)
          ! screen with oscillator strength
@@ -219,7 +220,7 @@ SUBROUTINE tddfpt_print_forces(qs_env, evects, evals, ostrength, print_section,
                ELSE
                   iw = -1
                END IF
-               CALL response_equation(qs_env, p_env, ex_env%cpmos, iw)
+               CALL response_equation(qs_env, p_env, ex_env%cpmos, iw, tdlr_section)
                !
                CALL get_qs_env(qs_env, dft_control=dft_control)
                IF (dft_control%qs_control%semi_empirical) THEN

diff --git a/tests/QS/regtest-stda-force/h2o_f01.inp b/tests/QS/regtest-stda-force/h2o_f01.inp
@@ -6,6 +6,12 @@
        &sTDA
           FRACTION 0.50
        &END sTDA
+       &LINRES
+         PRECONDITIONER FULL_SINGLE_INVERSE
+         MAX_ITER 20
+         EPS 1.0E-10
+         EPS_FILTER 1.0E-15
+       &END LINRES
        NSTATES      5
        MAX_ITER    50
        CONVERGENCE [eV] 1.0e-7

diff --git a/tests/QS/regtest-stda-force/h2o_f02.inp b/tests/QS/regtest-stda-force/h2o_f02.inp
@@ -6,6 +6,12 @@
        &sTDA
           FRACTION 0.50
        &END sTDA
+       &LINRES
+         PRECONDITIONER FULL_SINGLE_INVERSE
+         MAX_ITER 20
+         EPS 1.0E-10
+         EPS_FILTER 1.0E-15
+       &END LINRES
        NSTATES      5
        MAX_ITER    50
        CONVERGENCE [eV] 1.0e-7

diff --git a/tests/QS/regtest-stda/CH2O_s1.inp b/tests/QS/regtest-stda/CH2O_s1.inp
@@ -12,6 +12,12 @@
           FRACTION 0.0
           DO_EXCHANGE F
        &END sTDA
+       &LINRES
+         PRECONDITIONER FULL_SINGLE_INVERSE
+         MAX_ITER 20
+         EPS 1.0E-10
+         EPS_FILTER 1.0E-15
+       &END LINRES
        NSTATES      5
        MAX_ITER    50
        CONVERGENCE [eV] 1.0e-8

diff --git a/tests/QS/regtest-tddfpt-admm/h2o_a01.inp b/tests/QS/regtest-tddfpt-admm/h2o_a01.inp
@@ -9,6 +9,12 @@
          &XC_FUNCTIONAL PBE
          &END XC_FUNCTIONAL
        &END XC
+       &LINRES
+         PRECONDITIONER FULL_SINGLE_INVERSE
+         MAX_ITER 20
+         EPS 1.0E-10
+         EPS_FILTER 1.0E-15
+       &END LINRES
        NSTATES      2
        MAX_ITER    50
        CONVERGENCE [eV] 1.0e-6

diff --git a/tests/QS/regtest-tddfpt-admm/h2o_a02.inp b/tests/QS/regtest-tddfpt-admm/h2o_a02.inp
@@ -8,6 +8,12 @@
        NSTATES      2
        MAX_ITER    50
        CONVERGENCE [eV] 1.0e-5
+       &LINRES
+         PRECONDITIONER FULL_SINGLE_INVERSE
+         MAX_ITER 20
+         EPS 1.0E-10
+         EPS_FILTER 1.0E-15
+       &END LINRES
        RKS_TRIPLETS F
     &END TDDFPT
   &END PROPERTIES

diff --git a/tests/QS/regtest-tddfpt-force/h2o_f01.inp b/tests/QS/regtest-tddfpt-force/h2o_f01.inp
@@ -7,6 +7,12 @@
          &XC_FUNCTIONAL NONE
          &END XC_FUNCTIONAL
        &END XC
+       &LINRES
+         PRECONDITIONER FULL_SINGLE_INVERSE
+         MAX_ITER 20
+         EPS 1.0E-10
+         EPS_FILTER 1.0E-15
+       &END LINRES
        NSTATES      5
        MAX_ITER    50
        CONVERGENCE [eV] 1.0e-7

diff --git a/tests/QS/regtest-tddfpt-force/h2o_f02.inp b/tests/QS/regtest-tddfpt-force/h2o_f02.inp
@@ -7,6 +7,12 @@
          &XC_FUNCTIONAL PADE
          &END XC_FUNCTIONAL
        &END XC
+       &LINRES
+         PRECONDITIONER FULL_SINGLE_INVERSE
+         MAX_ITER 20
+         EPS 1.0E-10
+         EPS_FILTER 1.0E-15
+       &END LINRES
        NSTATES      5
        MAX_ITER    50
        CONVERGENCE [eV] 1.0e-7

diff --git a/tests/QS/regtest-tddfpt-force/h2o_f03.inp b/tests/QS/regtest-tddfpt-force/h2o_f03.inp
@@ -7,6 +7,12 @@
          &XC_FUNCTIONAL PBE
          &END XC_FUNCTIONAL
        &END XC
+       &LINRES
+         PRECONDITIONER FULL_SINGLE_INVERSE
+         MAX_ITER 20
+         EPS 1.0E-10
+         EPS_FILTER 1.0E-15
+       &END LINRES
        NSTATES      5
        MAX_ITER    50
        CONVERGENCE [eV] 1.0e-7

diff --git a/tests/QS/regtest-tddfpt-force/h2o_statefollow.inp b/tests/QS/regtest-tddfpt-force/h2o_statefollow.inp
@@ -7,6 +7,12 @@
        MAX_ITER    50
        CONVERGENCE [eV] 1.0e-7
        RKS_TRIPLETS F
+       &LINRES
+         PRECONDITIONER FULL_SINGLE_INVERSE
+         MAX_ITER 20
+         EPS 1.0E-10
+         EPS_FILTER 1.0E-15
+       &END LINRES
     &END TDDFPT
   &END PROPERTIES
   &DFT

diff --git a/tests/QS/regtest-tddfpt-prop/h2o_f01.inp b/tests/QS/regtest-tddfpt-prop/h2o_f01.inp
@@ -6,6 +6,10 @@
        NSTATES      5
        MAX_ITER    50
        CONVERGENCE [eV] 1.0e-7
+       &LINRES
+         PRECONDITIONER FULL_ALL
+         MAX_ITER 20
+       &END LINRES
        &PRINT
           &FORCES
              THRESHOLD 0.001