1!--------------------------------------------------------------------------------------------------! 2! CP2K: A general program to perform molecular dynamics simulations ! 3! Copyright (C) 2000 - 2019 CP2K developers group ! 4!--------------------------------------------------------------------------------------------------! 5 6! ************************************************************************************************** 7!> \brief Auxiliary routines for performing a constrained DFT SCF run with Quickstep. 8!> \par History 9!> - Separated some routines from qs_scf (03.2018) [Nico Holmberg] 10!> \author Nico Holmberg (03.2018) 11! ************************************************************************************************** 12MODULE qs_cdft_scf_utils 13 USE cp_control_types, ONLY: dft_control_type 14 USE cp_files, ONLY: close_file,& 15 get_unit_number,& 16 open_file 17 USE cp_log_handling, ONLY: cp_logger_create,& 18 cp_logger_type,& 19 cp_to_string 20 USE cp_para_types, ONLY: cp_para_env_type 21 USE input_constants, ONLY: & 22 broyden_type_1, broyden_type_1_explicit, broyden_type_1_explicit_ls, broyden_type_1_ls, & 23 broyden_type_2, broyden_type_2_explicit, broyden_type_2_explicit_ls, broyden_type_2_ls, & 24 jacobian_fd1, jacobian_fd1_backward, jacobian_fd1_central, jacobian_fd2, & 25 jacobian_fd2_backward, outer_scf_optimizer_broyden, outer_scf_optimizer_newton, & 26 outer_scf_optimizer_newton_ls 27 USE kinds, ONLY: default_path_length,& 28 dp 29 USE mathlib, ONLY: invert_matrix 30 USE qs_environment_types, ONLY: get_qs_env,& 31 qs_environment_type 32 USE qs_scf_types, ONLY: qs_scf_env_type 33 USE scf_control_types, ONLY: scf_control_type 34#include "./base/base_uses.f90" 35 36 IMPLICIT NONE 37 38 PRIVATE 39 40 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_cdft_scf_utils' 41 42 PUBLIC :: prepare_jacobian_stencil, build_diagonal_jacobian, & 43 restart_inverse_jacobian, print_inverse_jacobian, & 44 create_tmp_logger, initialize_inverse_jacobian 45 46CONTAINS 47 48! ************************************************************************************************** 49!> \brief Prepares the finite difference stencil for computing the Jacobian. The constraints 50!> are re-evaluated by perturbing each constraint. 51!> \param qs_env the qs_env where to build the Jacobian 52!> \param output_unit the output unit number 53!> \param nwork the number of perturbations to take in the negative direction 54!> \param pwork the number of perturbations to take in the positive direction 55!> \param coeff list of coefficients that determine how to sum up the various perturbations 56!> \param step_multiplier list of values that determine how large steps to take for each perturbatio 57!> \param dh total length of the interval to use for computing the finite difference derivatives 58!> \par History 59!> 03.2018 created [Nico Holmberg] 60! ************************************************************************************************** 61 SUBROUTINE prepare_jacobian_stencil(qs_env, output_unit, nwork, pwork, & 62 coeff, step_multiplier, dh) 63 TYPE(qs_environment_type), POINTER :: qs_env 64 INTEGER :: output_unit, nwork, pwork 65 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: coeff, step_multiplier, dh 66 67 CHARACTER(LEN=*), PARAMETER :: routineN = 'prepare_jacobian_stencil', & 68 routineP = moduleN//':'//routineN 69 70 CHARACTER(len=15) :: fmt_code 71 INTEGER :: ivar 72 TYPE(dft_control_type), POINTER :: dft_control 73 TYPE(qs_scf_env_type), POINTER :: scf_env 74 TYPE(scf_control_type), POINTER :: scf_control 75 76 NULLIFY (scf_env, scf_control, dft_control) 77 78 CPASSERT(ASSOCIATED(qs_env)) 79 CALL get_qs_env(qs_env, scf_env=scf_env, & 80 scf_control=scf_control, & 81 dft_control=dft_control) 82 83 IF (SIZE(scf_control%outer_scf%cdft_opt_control%jacobian_step) /= 1 .AND. & 84 SIZE(scf_control%outer_scf%cdft_opt_control%jacobian_step) /= SIZE(scf_env%outer_scf%variables, 1)) & 85 CALL cp_abort(__LOCATION__, & 86 cp_to_string(SIZE(scf_control%outer_scf%cdft_opt_control%jacobian_step))// & 87 " values passed to keyword JACOBIAN_STEP, expected 1 or "// & 88 cp_to_string(SIZE(scf_env%outer_scf%variables, 1))) 89 90 ALLOCATE (dh(SIZE(scf_env%outer_scf%variables, 1))) 91 IF (SIZE(dh) /= SIZE(scf_control%outer_scf%cdft_opt_control%jacobian_step)) THEN 92 DO ivar = 1, SIZE(dh) 93 dh(ivar) = scf_control%outer_scf%cdft_opt_control%jacobian_step(1) 94 END DO 95 ELSE 96 dh(:) = scf_control%outer_scf%cdft_opt_control%jacobian_step 97 END IF 98 99 SELECT CASE (scf_control%outer_scf%cdft_opt_control%jacobian_type) 100 CASE DEFAULT 101 CALL cp_abort(__LOCATION__, & 102 "Unknown Jacobian type: "// & 103 cp_to_string(scf_control%outer_scf%cdft_opt_control%jacobian_type)) 104 CASE (jacobian_fd1) 105 ! f'(x0) = [ f(x0+h) - f(x0) ] / h 106 nwork = 0 107 pwork = 1 108 ALLOCATE (coeff(nwork:pwork), step_multiplier(nwork:pwork)) 109 coeff(nwork) = -1.0_dp 110 coeff(pwork) = 1.0_dp 111 step_multiplier = 1.0_dp 112 CASE (jacobian_fd1_backward) 113 ! f'(x0) = [ f(x0) - f(x0-h) ] / h 114 nwork = -1 115 pwork = 0 116 ALLOCATE (coeff(nwork:pwork), step_multiplier(nwork:pwork)) 117 coeff(nwork) = -1.0_dp 118 coeff(pwork) = 1.0_dp 119 step_multiplier = -1.0_dp 120 CASE (jacobian_fd2) 121 ! f'(x0) = [ -f(x0+2h) + 4f(x0+h) - 3f(x0) ] / 2h 122 nwork = 0 123 pwork = 2 124 ALLOCATE (coeff(nwork:pwork), step_multiplier(nwork:pwork)) 125 coeff(0) = -3.0_dp 126 coeff(1) = 4.0_dp 127 coeff(2) = -1.0_dp 128 step_multiplier(0) = 0.0_dp 129 step_multiplier(1) = 1.0_dp 130 step_multiplier(2) = 2.0_dp 131 dh(:) = 2.0_dp*dh(:) 132 CASE (jacobian_fd2_backward) 133 ! f'(x0) = [ 3f(x0) - 4f(x0-h) + f(x0-2h) ] / 2h 134 nwork = -2 135 pwork = 0 136 ALLOCATE (coeff(nwork:pwork), step_multiplier(nwork:pwork)) 137 coeff(0) = 3.0_dp 138 coeff(-1) = -4.0_dp 139 coeff(-2) = 1.0_dp 140 step_multiplier(0) = 0.0_dp 141 step_multiplier(-1) = -1.0_dp 142 step_multiplier(-2) = -2.0_dp 143 dh(:) = 2.0_dp*dh(:) 144 CASE (jacobian_fd1_central) 145 ! f'(x0) = [ f(x0+h) - f(x0-h) ] / 2h 146 nwork = -1 147 pwork = 1 148 ALLOCATE (coeff(nwork:pwork), step_multiplier(nwork:pwork)) 149 coeff(0) = 0.0_dp 150 coeff(nwork) = -1.0_dp 151 coeff(pwork) = 1.0_dp 152 step_multiplier(0) = 0.0_dp 153 step_multiplier(nwork) = -1.0_dp 154 step_multiplier(pwork) = 1.0_dp 155 dh(:) = 2.0_dp*dh(:) 156 END SELECT 157 ! Print some info 158 IF (output_unit > 0) THEN 159 WRITE (output_unit, FMT="(/,A)") & 160 " ================================== JACOBIAN CALCULATION =================================" 161 WRITE (output_unit, FMT="(A)") & 162 " Evaluating inverse Jacobian using finite differences" 163 WRITE (output_unit, '(A,I10,A,I10)') & 164 " Energy evaluation: ", dft_control%qs_control%cdft_control%ienergy, & 165 ", CDFT SCF iteration: ", scf_env%outer_scf%iter_count 166 SELECT CASE (scf_control%outer_scf%cdft_opt_control%jacobian_type) 167 CASE (jacobian_fd1) 168 WRITE (output_unit, '(A)') " Type : First order forward difference" 169 CASE (jacobian_fd1_backward) 170 WRITE (output_unit, '(A)') " Type : First order backward difference" 171 CASE (jacobian_fd2) 172 WRITE (output_unit, '(A)') " Type : Second order forward difference" 173 CASE (jacobian_fd2_backward) 174 WRITE (output_unit, '(A)') " Type : Second order backward difference" 175 CASE (jacobian_fd1_central) 176 WRITE (output_unit, '(A)') " Type : First order central difference" 177 CASE DEFAULT 178 CALL cp_abort(__LOCATION__, "Unknown Jacobian type: "// & 179 cp_to_string(scf_control%outer_scf%cdft_opt_control%jacobian_type)) 180 END SELECT 181 IF (SIZE(scf_control%outer_scf%cdft_opt_control%jacobian_step) == 1) THEN 182 WRITE (output_unit, '(A,ES12.4)') " Step size : ", scf_control%outer_scf%cdft_opt_control%jacobian_step 183 ELSE 184 WRITE (output_unit, '(A,ES12.4,A)') & 185 " Step sizes : ", scf_control%outer_scf%cdft_opt_control%jacobian_step(1), ' (constraint 1)' 186 IF (SIZE(scf_control%outer_scf%cdft_opt_control%jacobian_step) < 10) THEN 187 fmt_code = '(ES34.4,A,I2,A)' 188 ELSE 189 fmt_code = '(ES34.4,A,I3,A)' 190 END IF 191 DO ivar = 2, SIZE(scf_control%outer_scf%cdft_opt_control%jacobian_step) 192 WRITE (output_unit, fmt_code) scf_control%outer_scf%cdft_opt_control%jacobian_step(ivar), " (constraint", ivar, ")" 193 END DO 194 END IF 195 END IF 196 197 END SUBROUTINE prepare_jacobian_stencil 198! ************************************************************************************************** 199!> \brief Builds a strictly diagonal inverse Jacobian from MD/SCF history. 200!> \param qs_env the qs_environment_type where to compute the Jacobian 201!> \param used_history flag that determines if history was actually used to prepare the Jacobian 202!> \par History 203!> 03.2018 created [Nico Holmberg] 204! ************************************************************************************************** 205 SUBROUTINE build_diagonal_jacobian(qs_env, used_history) 206 TYPE(qs_environment_type), POINTER :: qs_env 207 LOGICAL :: used_history 208 209 CHARACTER(LEN=*), PARAMETER :: routineN = 'build_diagonal_jacobian', & 210 routineP = moduleN//':'//routineN 211 212 INTEGER :: i, ihistory, nvar, outer_scf_ihistory 213 LOGICAL :: use_md_history 214 REAL(KIND=dp) :: inv_error 215 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: jacobian 216 REAL(KIND=dp), DIMENSION(:, :), POINTER :: gradient_history, inv_jacobian, & 217 variable_history 218 TYPE(qs_scf_env_type), POINTER :: scf_env 219 TYPE(scf_control_type), POINTER :: scf_control 220 221 NULLIFY (scf_control, scf_env) 222 223 CALL get_qs_env(qs_env, scf_env=scf_env, & 224 scf_control=scf_control, & 225 outer_scf_ihistory=outer_scf_ihistory) 226 ihistory = scf_env%outer_scf%iter_count 227 228 IF (outer_scf_ihistory .GE. 3 .AND. .NOT. used_history) THEN 229 ! First, lets try using the history from previous energy evaluations 230 CALL get_qs_env(qs_env, gradient_history=gradient_history, & 231 variable_history=variable_history) 232 nvar = SIZE(scf_env%outer_scf%variables, 1) 233 use_md_history = .TRUE. 234 ! Check that none of the history values are identical in which case we should try something different 235 DO i = 1, nvar 236 IF (ABS(variable_history(i, 2) - variable_history(i, 1)) .LT. 1.0E-12_dp) & 237 use_md_history = .FALSE. 238 END DO 239 IF (use_md_history) THEN 240 ALLOCATE (jacobian(nvar, nvar)) 241 DO i = 1, nvar 242 jacobian(i, i) = (gradient_history(i, 2) - gradient_history(i, 1))/ & 243 (variable_history(i, 2) - variable_history(i, 1)) 244 END DO 245 IF (.NOT. ASSOCIATED(scf_env%outer_scf%inv_jacobian)) & 246 ALLOCATE (scf_env%outer_scf%inv_jacobian(nvar, nvar)) 247 inv_jacobian => scf_env%outer_scf%inv_jacobian 248 CALL invert_matrix(jacobian, inv_jacobian, inv_error) 249 DEALLOCATE (jacobian) 250 ! Mark that an inverse Jacobian was just built and the next outer_loop_optimize should not perform 251 ! a Broyden update of it 252 scf_control%outer_scf%cdft_opt_control%broyden_update = .FALSE. 253 ! Mark that the MD history has been used and should not be reused anymore on this energy evaluation 254 used_history = .TRUE. 255 END IF 256 END IF 257 IF (ihistory .GE. 2 .AND. .NOT. ASSOCIATED(scf_env%outer_scf%inv_jacobian)) THEN 258 ! Next, try history from current SCF procedure 259 nvar = SIZE(scf_env%outer_scf%variables, 1) 260 IF (SIZE(scf_env%outer_scf%gradient, 2) .LT. 3) & 261 CALL cp_abort(__LOCATION__, & 262 "Keyword EXTRAPOLATION_ORDER in section OUTER_SCF must be greater than or equal "// & 263 "to 3 for optimizers that build the Jacobian from SCF history.") 264 ALLOCATE (jacobian(nvar, nvar)) 265 DO i = 1, nvar 266 jacobian(i, i) = (scf_env%outer_scf%gradient(i, ihistory) - scf_env%outer_scf%gradient(i, ihistory - 1))/ & 267 (scf_env%outer_scf%variables(i, ihistory) - scf_env%outer_scf%variables(i, ihistory - 1)) 268 END DO 269 IF (.NOT. ASSOCIATED(scf_env%outer_scf%inv_jacobian)) & 270 ALLOCATE (scf_env%outer_scf%inv_jacobian(nvar, nvar)) 271 inv_jacobian => scf_env%outer_scf%inv_jacobian 272 CALL invert_matrix(jacobian, inv_jacobian, inv_error) 273 DEALLOCATE (jacobian) 274 scf_control%outer_scf%cdft_opt_control%broyden_update = .FALSE. 275 ELSE 276 ! No history => will fall back to SD optimizer in outer_loop_optimize 277 END IF 278 279 END SUBROUTINE build_diagonal_jacobian 280 281! ************************************************************************************************** 282!> \brief Restarts the finite difference inverse Jacobian. 283!> \param qs_env the qs_environment_type where to compute the Jacobian 284!> \par History 285!> 03.2018 created [Nico Holmberg] 286! ************************************************************************************************** 287 SUBROUTINE restart_inverse_jacobian(qs_env) 288 TYPE(qs_environment_type), POINTER :: qs_env 289 290 CHARACTER(LEN=*), PARAMETER :: routineN = 'restart_inverse_jacobian', & 291 routineP = moduleN//':'//routineN 292 293 INTEGER :: i, iwork, j, nvar 294 REAL(KIND=dp), DIMENSION(:, :), POINTER :: inv_jacobian 295 TYPE(qs_scf_env_type), POINTER :: scf_env 296 TYPE(scf_control_type), POINTER :: scf_control 297 298 NULLIFY (scf_env, scf_control) 299 CPASSERT(ASSOCIATED(qs_env)) 300 CALL get_qs_env(qs_env, scf_env=scf_env, scf_control=scf_control) 301 302 CPASSERT(ASSOCIATED(scf_control%outer_scf%cdft_opt_control%jacobian_vector)) 303 nvar = SIZE(scf_env%outer_scf%variables, 1) 304 IF (SIZE(scf_control%outer_scf%cdft_opt_control%jacobian_vector) /= nvar**2) & 305 CALL cp_abort(__LOCATION__, & 306 "Too many or too few values defined for restarting inverse Jacobian.") 307 IF (.NOT. ASSOCIATED(scf_env%outer_scf%inv_jacobian)) & 308 ALLOCATE (scf_env%outer_scf%inv_jacobian(nvar, nvar)) 309 inv_jacobian => scf_env%outer_scf%inv_jacobian 310 iwork = 1 311 DO i = 1, nvar 312 DO j = 1, nvar 313 inv_jacobian(i, j) = scf_control%outer_scf%cdft_opt_control%jacobian_vector(iwork) 314 iwork = iwork + 1 315 END DO 316 END DO 317 DEALLOCATE (scf_control%outer_scf%cdft_opt_control%jacobian_vector) 318 scf_control%outer_scf%cdft_opt_control%jacobian_restart = .FALSE. 319 scf_control%outer_scf%cdft_opt_control%broyden_update = .FALSE. 320 scf_env%outer_scf%deallocate_jacobian = .FALSE. 321 322 END SUBROUTINE restart_inverse_jacobian 323 324! ************************************************************************************************** 325!> \brief Prints the finite difference inverse Jacobian to file 326!> \param logger the default IO logger 327!> \param inv_jacobian the inverse Jacobian matrix 328!> \param iter_count the iteration number 329!> \par History 330!> 03.2018 created [Nico Holmberg] 331! ************************************************************************************************** 332 SUBROUTINE print_inverse_jacobian(logger, inv_jacobian, iter_count) 333 TYPE(cp_logger_type), POINTER :: logger 334 REAL(KIND=dp), DIMENSION(:, :), INTENT(IN), & 335 POINTER :: inv_jacobian 336 INTEGER :: iter_count 337 338 CHARACTER(LEN=*), PARAMETER :: routineN = 'print_inverse_jacobian', & 339 routineP = moduleN//':'//routineN 340 341 CHARACTER(len=default_path_length) :: project_name 342 INTEGER :: i, j, lp, nvar, output_unit 343 344 nvar = SIZE(inv_jacobian, 1) 345 output_unit = get_unit_number() 346 project_name = logger%iter_info%project_name 347 lp = LEN_TRIM(project_name) 348 project_name(lp + 1:LEN(project_name)) = ".inverseJacobian" 349 CALL open_file(file_name=project_name, file_status="UNKNOWN", & 350 file_action="WRITE", file_position="APPEND", & 351 unit_number=output_unit) 352 WRITE (output_unit, FMT="(/,A)") "Inverse Jacobian matrix in row major order" 353 WRITE (output_unit, FMT="(A,I10)") "Iteration: ", iter_count 354 DO i = 1, nvar 355 DO j = 1, nvar 356 WRITE (output_unit, *) inv_jacobian(i, j) 357 END DO 358 END DO 359 CALL close_file(unit_number=output_unit) 360 361 END SUBROUTINE print_inverse_jacobian 362 363! ************************************************************************************************** 364!> \brief Creates a temporary logger for redirecting output to a new file 365!> \param para_env the para_env 366!> \param project_name the project basename 367!> \param suffix the suffix 368!> \param output_unit the default unit number for the newly created temporary logger 369!> \param tmp_logger pointer to the newly created temporary logger 370!> \par History 371!> 03.2018 created [Nico Holmberg] 372! ************************************************************************************************** 373 SUBROUTINE create_tmp_logger(para_env, project_name, suffix, output_unit, tmp_logger) 374 TYPE(cp_para_env_type), POINTER :: para_env 375 CHARACTER(len=*) :: project_name, suffix 376 INTEGER, INTENT(OUT) :: output_unit 377 TYPE(cp_logger_type), INTENT(OUT), POINTER :: tmp_logger 378 379 CHARACTER(LEN=*), PARAMETER :: routineN = 'create_tmp_logger', & 380 routineP = moduleN//':'//routineN 381 382 INTEGER :: lp 383 384 IF (para_env%ionode) THEN 385 lp = LEN_TRIM(project_name) 386 project_name(lp + 1:LEN(project_name)) = suffix 387 CALL open_file(file_name=project_name, file_status="UNKNOWN", & 388 file_action="WRITE", file_position="APPEND", & 389 unit_number=output_unit) 390 ELSE 391 output_unit = -1 392 END IF 393 CALL cp_logger_create(tmp_logger, & 394 para_env=para_env, & 395 default_global_unit_nr=output_unit, & 396 close_global_unit_on_dealloc=.FALSE.) 397 398 END SUBROUTINE create_tmp_logger 399 400! ************************************************************************************************** 401!> \brief Checks if the inverse Jacobian should be calculated and initializes the calculation 402!> \param scf_control the scf_control that holds the Jacobian settings 403!> \param scf_env the scf_env that holds the CDFT iteration information 404!> \param explicit_jacobian flag that determines if the finite difference Jacobian is needed 405!> \param should_build flag that determines if the Jacobian should be built 406!> \param used_history flag that determines if SCF history has been used to build a Jacobian 407!> \par History 408!> 03.2018 created [Nico Holmberg] 409! ************************************************************************************************** 410 SUBROUTINE initialize_inverse_jacobian(scf_control, scf_env, explicit_jacobian, & 411 should_build, used_history) 412 TYPE(scf_control_type), POINTER :: scf_control 413 TYPE(qs_scf_env_type), POINTER :: scf_env 414 LOGICAL :: explicit_jacobian, should_build, & 415 used_history 416 417 CHARACTER(LEN=*), PARAMETER :: routineN = 'initialize_inverse_jacobian', & 418 routineP = moduleN//':'//routineN 419 420 CPASSERT(ASSOCIATED(scf_control)) 421 CPASSERT(ASSOCIATED(scf_env)) 422 423 SELECT CASE (scf_control%outer_scf%optimizer) 424 CASE DEFAULT 425 CPABORT("Noncompatible optimizer requested.") 426 CASE (outer_scf_optimizer_newton, outer_scf_optimizer_newton_ls) 427 CPASSERT(ASSOCIATED(scf_control%outer_scf%cdft_opt_control)) 428 scf_control%outer_scf%cdft_opt_control%build_jacobian = .TRUE. 429 explicit_jacobian = .TRUE. 430 CASE (outer_scf_optimizer_broyden) 431 CPASSERT(ASSOCIATED(scf_control%outer_scf%cdft_opt_control)) 432 SELECT CASE (scf_control%outer_scf%cdft_opt_control%broyden_type) 433 CASE (broyden_type_1, broyden_type_2, broyden_type_1_ls, broyden_type_2_ls) 434 scf_control%outer_scf%cdft_opt_control%build_jacobian = .TRUE. 435 explicit_jacobian = .FALSE. 436 CASE (broyden_type_1_explicit, broyden_type_2_explicit, broyden_type_1_explicit_ls, broyden_type_2_explicit_ls) 437 scf_control%outer_scf%cdft_opt_control%build_jacobian = .TRUE. 438 explicit_jacobian = .TRUE. 439 END SELECT 440 END SELECT 441 IF (scf_control%outer_scf%cdft_opt_control%build_jacobian) THEN 442 ! Reset counter 443 IF (scf_env%outer_scf%iter_count == 1) scf_control%outer_scf%cdft_opt_control%ijacobian(1) = 0 444 ! Check if an old Jacobian can be reused avoiding a rebuild 445 IF (ASSOCIATED(scf_env%outer_scf%inv_jacobian)) THEN 446 ! Rebuild if number of previous energy evaluations exceeds limit 447 IF (scf_control%outer_scf%cdft_opt_control%ijacobian(2) .GE. & 448 scf_control%outer_scf%cdft_opt_control%jacobian_freq(2) .AND. .NOT. used_history .AND. & 449 scf_control%outer_scf%cdft_opt_control%jacobian_freq(2) > 0) THEN 450 should_build = .TRUE. 451 ! Zero the corresponding counters 452 scf_control%outer_scf%cdft_opt_control%ijacobian(:) = 0 453 ! Rebuild if number of previous SCF iterations exceeds limit (on this energy eval) 454 ELSE IF (scf_control%outer_scf%cdft_opt_control%ijacobian(1) .GE. & 455 scf_control%outer_scf%cdft_opt_control%jacobian_freq(1) .AND. & 456 scf_control%outer_scf%cdft_opt_control%jacobian_freq(1) > 0) THEN 457 should_build = .TRUE. 458 ! Zero the corresponding counter 459 scf_control%outer_scf%cdft_opt_control%ijacobian(1) = 0 460 END IF 461 IF (should_build) DEALLOCATE (scf_env%outer_scf%inv_jacobian) 462 ELSE 463 should_build = .TRUE. 464 ! Zero the counter 465 scf_control%outer_scf%cdft_opt_control%ijacobian(:) = 0 466 END IF 467 END IF 468 469 END SUBROUTINE initialize_inverse_jacobian 470 471END MODULE qs_cdft_scf_utils 472