1!--------------------------------------------------------------------------------------------------! 2! CP2K: A general program to perform molecular dynamics simulations ! 3! Copyright (C) 2000 - 2019 CP2K developers group ! 4!--------------------------------------------------------------------------------------------------! 5 6! ************************************************************************************************** 7!> \brief Control parameters for optimizers that work with CDFT constraints 8!> \par History 9!> separated from scf_control_types [03.2018] 10!> \author Nico Holmberg [03.2018] 11! ************************************************************************************************** 12MODULE qs_cdft_opt_types 13 14 USE input_constants, ONLY: & 15 broyden_type_1, broyden_type_1_explicit, broyden_type_1_explicit_ls, broyden_type_1_ls, & 16 broyden_type_2, broyden_type_2_explicit, broyden_type_2_explicit_ls, broyden_type_2_ls, & 17 outer_scf_optimizer_broyden, outer_scf_optimizer_newton, outer_scf_optimizer_newton_ls 18 USE input_section_types, ONLY: section_vals_get_subs_vals,& 19 section_vals_type,& 20 section_vals_val_get 21 USE kinds, ONLY: dp 22#include "./base/base_uses.f90" 23 24 IMPLICIT NONE 25 26 PRIVATE 27 28 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_cdft_opt_types' 29 30 ! Public data types 31 32 PUBLIC :: cdft_opt_type 33 34 ! Public subroutines 35 36 PUBLIC :: cdft_opt_type_create, & 37 cdft_opt_type_release, & 38 cdft_opt_type_read, & 39 cdft_opt_type_write, & 40 cdft_opt_type_copy 41 42! ************************************************************************************************** 43!> \brief contains the parameters needed by CDFT specific optimizers 44!> \param build_jacobian logical which determines if the inverse Jacobian should be computed 45!> \param jacobian_step the step size for calculating the finite difference Jacobian 46!> \param newton_step the step size used by the Newton optimizer with values between 0 and 1 47!> \param newton_step_save permanent copy of the above 48!> \param jacobian_type the finite difference scheme to compute the Jacobian 49!> \param broyden_type the variant of Broyden's method to use 50!> \param jacobian_freq control parameters defining how often the Jacobian is built 51!> \param ijacobian counter to track how many SCF iterations/energy evaluations have passed since 52!> the last Jacobian rebuild 53!> \param broyden_update logical which determines if a Broyden update is needed 54!> \param max_ls the maximum number of backtracking line search steps to perform 55!> \param continue_ls continue line search until max steps are reached or until the gradient 56!> no longer decreases 57!> \param factor_ls line search parameter used in generating a new step size 58!> \par History 59!> created [03.2018] 60!> \author Nico Holmberg [03.2018] 61! ************************************************************************************************** 62 63 TYPE cdft_opt_type 64 LOGICAL :: build_jacobian 65 LOGICAL :: broyden_update 66 LOGICAL :: continue_ls 67 LOGICAL :: jacobian_restart 68 REAL(KIND=dp) :: newton_step 69 REAL(KIND=dp) :: newton_step_save 70 REAL(KIND=dp) :: factor_ls 71 REAL(KIND=dp), DIMENSION(:), & 72 ALLOCATABLE :: jacobian_step 73 REAL(KIND=dp), DIMENSION(:), & 74 POINTER :: jacobian_vector 75 INTEGER :: jacobian_type 76 INTEGER :: broyden_type 77 INTEGER :: jacobian_freq(2) 78 INTEGER :: ijacobian(2) 79 INTEGER :: max_ls 80 END TYPE cdft_opt_type 81 82CONTAINS 83 84! ************************************************************************************************** 85!> \brief allocates and initializes the CDFT optimizer control object with default values 86!> \param cdft_opt_control the object to initialize 87!> \par History 88!> 03.2018 created [Nico Holmberg] 89!> \author Nico Holmberg 90! ************************************************************************************************** 91 SUBROUTINE cdft_opt_type_create(cdft_opt_control) 92 93 TYPE(cdft_opt_type), POINTER :: cdft_opt_control 94 95 CHARACTER(LEN=*), PARAMETER :: routineN = 'cdft_opt_type_create', & 96 routineP = moduleN//':'//routineN 97 98 INTEGER :: handle 99 100 CALL timeset(routineN, handle) 101 102 CPASSERT(.NOT. ASSOCIATED(cdft_opt_control)) 103 ALLOCATE (cdft_opt_control) 104 105 ! Load the default values 106 107 cdft_opt_control%jacobian_type = -1 108 cdft_opt_control%broyden_type = -1 109 cdft_opt_control%jacobian_freq(:) = 1 110 cdft_opt_control%newton_step = 1.0_dp 111 cdft_opt_control%newton_step_save = 1.0_dp 112 cdft_opt_control%factor_ls = 0.5_dp 113 cdft_opt_control%ijacobian(:) = 0 114 cdft_opt_control%max_ls = 0 115 cdft_opt_control%build_jacobian = .FALSE. 116 cdft_opt_control%broyden_update = .FALSE. 117 cdft_opt_control%continue_ls = .FALSE. 118 cdft_opt_control%jacobian_restart = .FALSE. 119 NULLIFY (cdft_opt_control%jacobian_vector) 120 121 CALL timestop(handle) 122 123 END SUBROUTINE cdft_opt_type_create 124 125! ************************************************************************************************** 126!> \brief releases the CDFT optimizer control object 127!> \param cdft_opt_control the object to release 128!> \par History 129!> 03.2018 created [Nico Holmberg] 130!> \author Nico Holmberg 131! ************************************************************************************************** 132 SUBROUTINE cdft_opt_type_release(cdft_opt_control) 133 134 TYPE(cdft_opt_type), POINTER :: cdft_opt_control 135 136 CHARACTER(LEN=*), PARAMETER :: routineN = 'cdft_opt_type_release', & 137 routineP = moduleN//':'//routineN 138 139 IF (ASSOCIATED(cdft_opt_control)) THEN 140 IF (ASSOCIATED(cdft_opt_control%jacobian_vector)) & 141 DEALLOCATE (cdft_opt_control%jacobian_vector) 142 IF (ALLOCATED(cdft_opt_control%jacobian_step)) & 143 DEALLOCATE (cdft_opt_control%jacobian_step) 144 145 DEALLOCATE (cdft_opt_control) 146 END IF 147 148 NULLIFY (cdft_opt_control) 149 150 END SUBROUTINE cdft_opt_type_release 151 152! ************************************************************************************************** 153!> \brief reads the parameters of the CDFT optimizer type 154!> \param cdft_opt_control the object that will contain the values read 155!> \param inp_section the input section that contains the values that are read 156!> \par History 157!> 03.2018 created [Nico Holmberg] 158!> \author Nico Holmberg 159! ************************************************************************************************** 160 SUBROUTINE cdft_opt_type_read(cdft_opt_control, inp_section) 161 162 TYPE(cdft_opt_type), POINTER :: cdft_opt_control 163 TYPE(section_vals_type), POINTER :: inp_section 164 165 CHARACTER(LEN=*), PARAMETER :: routineN = 'cdft_opt_type_read', & 166 routineP = moduleN//':'//routineN 167 168 INTEGER :: handle 169 INTEGER, DIMENSION(:), POINTER :: tmplist 170 LOGICAL :: exists 171 REAL(KIND=dp), DIMENSION(:), POINTER :: rtmplist 172 TYPE(section_vals_type), POINTER :: cdft_opt_section 173 174 CALL timeset(routineN, handle) 175 176 CPASSERT(ASSOCIATED(cdft_opt_control)) 177 cdft_opt_section => section_vals_get_subs_vals(inp_section, "CDFT_OPT") 178 179 CALL section_vals_val_get(cdft_opt_section, "MAX_LS", & 180 i_val=cdft_opt_control%max_ls) 181 CALL section_vals_val_get(cdft_opt_section, "JACOBIAN_TYPE", & 182 i_val=cdft_opt_control%jacobian_type) 183 CALL section_vals_val_get(cdft_opt_section, "JACOBIAN_STEP", & 184 r_vals=rtmplist) 185 ALLOCATE (cdft_opt_control%jacobian_step(SIZE(rtmplist))) 186 cdft_opt_control%jacobian_step(:) = rtmplist 187 CALL section_vals_val_get(cdft_opt_section, "BROYDEN_TYPE", & 188 i_val=cdft_opt_control%broyden_type) 189 CALL section_vals_val_get(cdft_opt_section, "CONTINUE_LS", & 190 l_val=cdft_opt_control%continue_ls) 191 CALL section_vals_val_get(cdft_opt_section, "FACTOR_LS", & 192 r_val=cdft_opt_control%factor_ls) 193 IF (cdft_opt_control%factor_ls .LE. 0.0_dp .OR. & 194 cdft_opt_control%factor_ls .GE. 1.0_dp) & 195 CALL cp_abort(__LOCATION__, & 196 "Keyword FACTOR_LS must be between 0.0 and 1.0.") 197 CALL section_vals_val_get(cdft_opt_section, "JACOBIAN_FREQ", explicit=exists) 198 IF (exists) THEN 199 CALL section_vals_val_get(cdft_opt_section, "JACOBIAN_FREQ", & 200 i_vals=tmplist) 201 IF (SIZE(tmplist) /= 2) & 202 CALL cp_abort(__LOCATION__, & 203 "Keyword JACOBIAN_FREQ takes exactly two input values.") 204 IF (ANY(tmplist .LT. 0)) & 205 CALL cp_abort(__LOCATION__, & 206 "Keyword JACOBIAN_FREQ takes only positive values.") 207 IF (ALL(tmplist .EQ. 0)) & 208 CALL cp_abort(__LOCATION__, & 209 "Both values to keyword JACOBIAN_FREQ cannot be zero.") 210 cdft_opt_control%jacobian_freq(:) = tmplist(1:2) 211 END IF 212 CALL section_vals_val_get(cdft_opt_section, "JACOBIAN_RESTART", & 213 l_val=cdft_opt_control%jacobian_restart) 214 IF (cdft_opt_control%jacobian_restart) THEN 215 CALL section_vals_val_get(cdft_opt_section, "JACOBIAN_VECTOR", & 216 r_vals=rtmplist) 217 ALLOCATE (cdft_opt_control%jacobian_vector(SIZE(rtmplist))) 218 cdft_opt_control%jacobian_vector = rtmplist 219 END IF 220 221 CALL timestop(handle) 222 223 END SUBROUTINE cdft_opt_type_read 224 225! ************************************************************************************************** 226!> \brief writes information about the CDFT optimizer object 227!> \param cdft_opt_control the CDFT optimizer object 228!> \param optimizer the type of optimizer to use 229!> \param output_unit the output unit handle 230!> \par History 231!> 03.2018 created [Nico Holmberg] 232!> \author Nico Holmberg 233! ************************************************************************************************** 234 SUBROUTINE cdft_opt_type_write(cdft_opt_control, optimizer, output_unit) 235 TYPE(cdft_opt_type), POINTER :: cdft_opt_control 236 INTEGER :: optimizer, output_unit 237 238 CHARACTER(LEN=*), PARAMETER :: routineN = 'cdft_opt_type_write', & 239 routineP = moduleN//':'//routineN 240 241 CPASSERT(ASSOCIATED(cdft_opt_control)) 242 243 SELECT CASE (optimizer) 244 CASE DEFAULT 245 ! Do nothing 246 CASE (outer_scf_optimizer_broyden) 247 WRITE (output_unit, '(T3,A)') "Optimization with Broyden's method" 248 SELECT CASE (cdft_opt_control%broyden_type) 249 CASE (broyden_type_1) 250 WRITE (output_unit, '(A)') " variant : 1st method" 251 CASE (broyden_type_1_explicit) 252 WRITE (output_unit, '(A)') " variant : 1st method with explicit initial Jacobian" 253 CASE (broyden_type_1_ls) 254 WRITE (output_unit, '(A)') " variant : 1st method with backtracking line search" 255 CASE (broyden_type_1_explicit_ls) 256 WRITE (output_unit, '(A)') & 257 " variant : 1st method with explicit initial Jacobian" 258 WRITE (output_unit, '(A)') & 259 " and backtracking line search" 260 CASE (broyden_type_2) 261 WRITE (output_unit, '(A)') " variant : 2nd method" 262 CASE (broyden_type_2_explicit) 263 WRITE (output_unit, '(A)') " variant : 2nd method with explicit initial Jacobian" 264 CASE (broyden_type_2_ls) 265 WRITE (output_unit, '(A)') " variant : 2nd method with backtracking line search" 266 CASE (broyden_type_2_explicit_ls) 267 WRITE (output_unit, '(A)') & 268 " variant : 2nd method with explicit initial Jacobian" 269 WRITE (output_unit, '(A)') & 270 " and backtracking line search" 271 END SELECT 272 CASE (outer_scf_optimizer_newton) 273 WRITE (output_unit, '(T3,A)') "Optimization with Newton's method" 274 CASE (outer_scf_optimizer_newton_ls) 275 WRITE (output_unit, '(T3,A)') "Optimization with Newton's method using backtracking line search" 276 END SELECT 277 SELECT CASE (optimizer) 278 CASE DEFAULT 279 ! Do nothing 280 CASE (outer_scf_optimizer_broyden, outer_scf_optimizer_newton, outer_scf_optimizer_newton_ls) 281 IF (cdft_opt_control%jacobian_freq(2) > 0) THEN 282 WRITE (output_unit, '(T6,A,I4,A)') & 283 "The Jacobian is restarted every ", cdft_opt_control%jacobian_freq(2), " energy evaluation" 284 IF (cdft_opt_control%jacobian_freq(1) > 0) & 285 WRITE (output_unit, '(T29,A,I4,A)') & 286 "or every ", cdft_opt_control%jacobian_freq(1), " CDFT SCF iteration" 287 ELSE 288 WRITE (output_unit, '(T6,A,I4,A)') & 289 "The Jacobian is restarted every ", cdft_opt_control%jacobian_freq(1), " CDFT SCF iteration" 290 END IF 291 WRITE (output_unit, '(T3,A,F8.4)') & 292 "Optimizer step size: ", cdft_opt_control%newton_step_save 293 END SELECT 294 295 END SUBROUTINE cdft_opt_type_write 296 297! ************************************************************************************************** 298!> \brief copies settings between two CDFT optimizer control objects retaining both 299!> \param new the object where to copy the settings 300!> \param old the object from where to copy the settings 301!> \par History 302!> 03.2018 created [Nico Holmberg] 303!> \author Nico Holmberg 304! ************************************************************************************************** 305 SUBROUTINE cdft_opt_type_copy(new, old) 306 307 TYPE(cdft_opt_type), POINTER :: new, old 308 309 CHARACTER(LEN=*), PARAMETER :: routineN = 'cdft_opt_type_copy', & 310 routineP = moduleN//':'//routineN 311 312 INTEGER :: handle 313 314 ! Do nothing if cdft_opt_type is not allocated 315 ! this happens if CDFT is performed with an optimizer other than Broyden/Newton 316 IF (.NOT. ASSOCIATED(old)) RETURN 317 318 CALL timeset(routineN, handle) 319 320 IF (.NOT. ASSOCIATED(new)) CALL cdft_opt_type_create(new) 321 new%max_ls = old%max_ls 322 new%continue_ls = old%continue_ls 323 new%factor_ls = old%factor_ls 324 new%jacobian_type = old%jacobian_type 325 new%jacobian_freq(:) = old%jacobian_freq(:) 326 new%newton_step = old%newton_step 327 new%newton_step_save = old%newton_step_save 328 new%ijacobian(:) = old%ijacobian(:) 329 new%build_jacobian = old%build_jacobian 330 new%broyden_type = old%broyden_type 331 new%broyden_update = old%broyden_update 332 IF (ALLOCATED(new%jacobian_step)) DEALLOCATE (new%jacobian_step) 333 ALLOCATE (new%jacobian_step(SIZE(old%jacobian_step))) 334 new%jacobian_step(:) = old%jacobian_step 335 IF (old%jacobian_restart) THEN 336 ! Transfer restart vector for inverse Jacobian matrix 337 ! (qs_calculate_inverse_jacobian handles deallocation of transferred vector) 338 new%jacobian_restart = .TRUE. 339 ALLOCATE (new%jacobian_vector(SIZE(old%jacobian_vector))) 340 new%jacobian_vector = old%jacobian_vector 341 DEALLOCATE (old%jacobian_vector) 342 old%jacobian_restart = .FALSE. 343 END IF 344 345 CALL timestop(handle) 346 347 END SUBROUTINE cdft_opt_type_copy 348 349END MODULE qs_cdft_opt_types 350