1!--------------------------------------------------------------------------------------------------! 2! CP2K: A general program to perform molecular dynamics simulations ! 3! Copyright (C) 2000 - 2020 CP2K developers group ! 4!--------------------------------------------------------------------------------------------------! 5 6! ************************************************************************************************** 7!> \brief Routines for all ALMO-based SCF methods 8!> 'RZK-warning' marks unresolved issues 9!> \par History 10!> 2011.05 created [Rustam Z Khaliullin] 11!> \author Rustam Z Khaliullin 12! ************************************************************************************************** 13MODULE almo_scf 14 USE almo_scf_methods, ONLY: almo_scf_p_blk_to_t_blk,& 15 almo_scf_t_rescaling,& 16 almo_scf_t_to_proj,& 17 distribute_domains,& 18 orthogonalize_mos 19 USE almo_scf_optimizer, ONLY: almo_scf_block_diagonal,& 20 almo_scf_construct_nlmos,& 21 almo_scf_xalmo_eigensolver,& 22 almo_scf_xalmo_pcg,& 23 almo_scf_xalmo_trustr 24 USE almo_scf_qs, ONLY: almo_dm_to_almo_ks,& 25 almo_scf_construct_quencher,& 26 calculate_w_matrix_almo,& 27 construct_qs_mos,& 28 init_almo_ks_matrix_via_qs,& 29 matrix_almo_create,& 30 matrix_qs_to_almo 31 USE almo_scf_types, ONLY: almo_mat_dim_aobasis,& 32 almo_mat_dim_occ,& 33 almo_mat_dim_virt,& 34 almo_mat_dim_virt_disc,& 35 almo_mat_dim_virt_full,& 36 almo_scf_env_type,& 37 optimizer_options_type,& 38 print_optimizer_options 39 USE atomic_kind_types, ONLY: atomic_kind_type 40 USE bibliography, ONLY: Khaliullin2013,& 41 Kolafa2004,& 42 Scheiber2018,& 43 Staub2019,& 44 cite_reference 45 USE cp_blacs_env, ONLY: cp_blacs_env_release,& 46 cp_blacs_env_retain 47 USE cp_control_types, ONLY: dft_control_type 48 USE cp_dbcsr_diag, ONLY: cp_dbcsr_syevd 49 USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm 50 USE cp_fm_types, ONLY: cp_fm_type 51 USE cp_log_handling, ONLY: cp_get_default_logger,& 52 cp_logger_get_default_unit_nr,& 53 cp_logger_type 54 USE cp_para_env, ONLY: cp_para_env_release,& 55 cp_para_env_retain 56 USE cp_para_types, ONLY: cp_para_env_type 57 USE dbcsr_api, ONLY: & 58 dbcsr_add, dbcsr_add_on_diag, dbcsr_binary_read, dbcsr_checksum, dbcsr_copy, dbcsr_create, & 59 dbcsr_distribution_get, dbcsr_distribution_type, dbcsr_filter, dbcsr_finalize, & 60 dbcsr_get_info, dbcsr_get_stored_coordinates, dbcsr_init_random, & 61 dbcsr_iterator_blocks_left, dbcsr_iterator_next_block, dbcsr_iterator_start, & 62 dbcsr_iterator_stop, dbcsr_iterator_type, dbcsr_multiply, dbcsr_nblkcols_total, & 63 dbcsr_nblkrows_total, dbcsr_p_type, dbcsr_release, dbcsr_reserve_block2d, dbcsr_scale, & 64 dbcsr_set, dbcsr_type, dbcsr_type_no_symmetry, dbcsr_type_symmetric, dbcsr_work_create 65 USE domain_submatrix_methods, ONLY: init_submatrices,& 66 release_submatrices 67 USE input_constants, ONLY: & 68 almo_deloc_none, almo_deloc_scf, almo_deloc_x, almo_deloc_x_then_scf, & 69 almo_deloc_xalmo_1diag, almo_deloc_xalmo_scf, almo_deloc_xalmo_x, almo_deloc_xk, & 70 almo_domain_layout_molecular, almo_mat_distr_atomic, almo_mat_distr_molecular, & 71 almo_scf_diag, almo_scf_dm_sign, almo_scf_pcg, almo_scf_skip, almo_scf_trustr, & 72 atomic_guess, molecular_guess, optimizer_diis, optimizer_lin_eq_pcg, optimizer_pcg, & 73 optimizer_trustr, restart_guess, smear_fermi_dirac, virt_full, virt_number, virt_occ_size, & 74 xalmo_case_block_diag, xalmo_case_fully_deloc, xalmo_case_normal, xalmo_trial_r0_out 75 USE input_section_types, ONLY: section_vals_type 76 USE iterate_matrix, ONLY: invert_Hotelling,& 77 matrix_sqrt_Newton_Schulz 78 USE kinds, ONLY: default_path_length,& 79 dp 80 USE mathlib, ONLY: binomial 81 USE message_passing, ONLY: mp_sum 82 USE molecule_types, ONLY: get_molecule_set_info,& 83 molecule_type 84 USE mscfg_types, ONLY: get_matrix_from_submatrices,& 85 molecular_scf_guess_env_type 86 USE particle_types, ONLY: particle_type 87 USE qs_environment_types, ONLY: get_qs_env,& 88 qs_environment_type 89 USE qs_initial_guess, ONLY: calculate_atomic_block_dm,& 90 calculate_mopac_dm 91 USE qs_kind_types, ONLY: qs_kind_type 92 USE qs_mo_types, ONLY: get_mo_set,& 93 mo_set_p_type 94 USE qs_rho_types, ONLY: qs_rho_get,& 95 qs_rho_type 96 USE qs_scf_post_scf, ONLY: qs_scf_compute_properties 97 USE qs_scf_types, ONLY: qs_scf_env_type,& 98 scf_env_release 99#include "./base/base_uses.f90" 100 101 IMPLICIT NONE 102 103 PRIVATE 104 105 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'almo_scf' 106 107 PUBLIC :: almo_entry_scf 108 109 LOGICAL, PARAMETER :: debug_mode = .FALSE. 110 LOGICAL, PARAMETER :: safe_mode = .FALSE. 111 112CONTAINS 113 114! ************************************************************************************************** 115!> \brief The entry point into ALMO SCF routines 116!> \param qs_env pointer to the QS environment 117!> \param calc_forces calculate forces? 118!> \par History 119!> 2011.05 created [Rustam Z Khaliullin] 120!> \author Rustam Z Khaliullin 121! ************************************************************************************************** 122 SUBROUTINE almo_entry_scf(qs_env, calc_forces) 123 TYPE(qs_environment_type), POINTER :: qs_env 124 LOGICAL, INTENT(IN) :: calc_forces 125 126 CHARACTER(len=*), PARAMETER :: routineN = 'almo_entry_scf', routineP = moduleN//':'//routineN 127 128 INTEGER :: handle 129 TYPE(almo_scf_env_type), POINTER :: almo_scf_env 130 131 CALL timeset(routineN, handle) 132 133 CALL cite_reference(Khaliullin2013) 134 135 ! get a pointer to the almo environment 136 CALL get_qs_env(qs_env, almo_scf_env=almo_scf_env) 137 138 ! initialize scf 139 CALL almo_scf_init(qs_env, almo_scf_env, calc_forces) 140 141 ! create the initial guess for ALMOs 142 CALL almo_scf_initial_guess(qs_env, almo_scf_env) 143 144 ! perform SCF for block diagonal ALMOs 145 CALL almo_scf_main(qs_env, almo_scf_env) 146 147 ! allow electron delocalization 148 CALL almo_scf_delocalization(qs_env, almo_scf_env) 149 150 ! construct NLMOs 151 CALL construct_nlmos(qs_env, almo_scf_env) 152 153 ! electron correlation methods 154 !CALL almo_correlation_main(qs_env,almo_scf_env) 155 156 ! do post scf processing 157 CALL almo_scf_post(qs_env, almo_scf_env) 158 159 ! clean up the mess 160 CALL almo_scf_clean_up(almo_scf_env) 161 162 CALL timestop(handle) 163 164 END SUBROUTINE almo_entry_scf 165 166! ************************************************************************************************** 167!> \brief Initialization of the almo_scf_env_type. 168!> \param qs_env ... 169!> \param almo_scf_env ... 170!> \param calc_forces ... 171!> \par History 172!> 2011.05 created [Rustam Z Khaliullin] 173!> 2018.09 smearing support [Ruben Staub] 174!> \author Rustam Z Khaliullin 175! ************************************************************************************************** 176 SUBROUTINE almo_scf_init(qs_env, almo_scf_env, calc_forces) 177 TYPE(qs_environment_type), POINTER :: qs_env 178 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 179 LOGICAL, INTENT(IN) :: calc_forces 180 181 CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_init', routineP = moduleN//':'//routineN 182 183 INTEGER :: ao, handle, i, iao, idomain, ispin, & 184 multip, naos, natoms, ndomains, nelec, & 185 nelec_a, nelec_b, nmols, nspins, & 186 unit_nr 187 TYPE(cp_logger_type), POINTER :: logger 188 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s 189 TYPE(dft_control_type), POINTER :: dft_control 190 TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set 191 TYPE(section_vals_type), POINTER :: input 192 193 CALL timeset(routineN, handle) 194 195 ! define the output_unit 196 logger => cp_get_default_logger() 197 IF (logger%para_env%ionode) THEN 198 unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.) 199 ELSE 200 unit_nr = -1 201 ENDIF 202 203 ! set optimizers' types 204 almo_scf_env%opt_block_diag_diis%optimizer_type = optimizer_diis 205 almo_scf_env%opt_block_diag_pcg%optimizer_type = optimizer_pcg 206 almo_scf_env%opt_xalmo_diis%optimizer_type = optimizer_diis 207 almo_scf_env%opt_xalmo_pcg%optimizer_type = optimizer_pcg 208 almo_scf_env%opt_xalmo_trustr%optimizer_type = optimizer_trustr 209 almo_scf_env%opt_nlmo_pcg%optimizer_type = optimizer_pcg 210 almo_scf_env%opt_xalmo_newton_pcg_solver%optimizer_type = optimizer_lin_eq_pcg 211 212 ! get info from the qs_env 213 CALL get_qs_env(qs_env, & 214 nelectron_total=almo_scf_env%nelectrons_total, & 215 matrix_s=matrix_s, & 216 dft_control=dft_control, & 217 molecule_set=molecule_set, & 218 input=input, & 219 has_unit_metric=almo_scf_env%orthogonal_basis, & 220 para_env=almo_scf_env%para_env, & 221 blacs_env=almo_scf_env%blacs_env, & 222 nelectron_spin=almo_scf_env%nelectrons_spin) 223 CALL cp_para_env_retain(almo_scf_env%para_env) 224 CALL cp_blacs_env_retain(almo_scf_env%blacs_env) 225 226 ! copy basic quantities 227 almo_scf_env%nspins = dft_control%nspins 228 almo_scf_env%natoms = dbcsr_nblkrows_total(matrix_s(1)%matrix) 229 almo_scf_env%nmolecules = SIZE(molecule_set) 230 CALL dbcsr_get_info(matrix_s(1)%matrix, nfullrows_total=naos) 231 almo_scf_env%naos = naos 232 233 !! retrieve smearing parameters, and check compatibility of methods requested 234 almo_scf_env%smear = dft_control%smear 235 IF (almo_scf_env%smear) THEN 236 CALL cite_reference(Staub2019) 237 IF ((almo_scf_env%almo_update_algorithm .NE. almo_scf_diag) .OR. & 238 ((almo_scf_env%deloc_method .NE. almo_deloc_none) .AND. & 239 (almo_scf_env%xalmo_update_algorithm .NE. almo_scf_diag))) THEN 240 CPABORT("ALMO smearing is currently implemented for DIAG algorithm only") 241 END IF 242 IF (qs_env%scf_control%smear%method .NE. smear_fermi_dirac) THEN 243 CPABORT("Only Fermi-Dirac smearing is currently compatible with ALMO") 244 END IF 245 almo_scf_env%smear_e_temp = qs_env%scf_control%smear%electronic_temperature 246 IF ((almo_scf_env%mat_distr_aos .NE. almo_mat_distr_molecular) .OR. & 247 (almo_scf_env%domain_layout_mos .NE. almo_domain_layout_molecular)) THEN 248 CPABORT("ALMO smearing was designed to work with molecular fragments only") 249 END IF 250 END IF 251 252 ! convenient local varibales 253 nspins = almo_scf_env%nspins 254 nmols = almo_scf_env%nmolecules 255 natoms = almo_scf_env%natoms 256 257 ! Define groups: either atomic or molecular 258 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular) THEN 259 almo_scf_env%ndomains = almo_scf_env%nmolecules 260 ELSE 261 almo_scf_env%ndomains = almo_scf_env%natoms 262 ENDIF 263 264 ! allocate domain descriptors 265 ndomains = almo_scf_env%ndomains 266 ALLOCATE (almo_scf_env%domain_index_of_atom(natoms)) 267 ALLOCATE (almo_scf_env%domain_index_of_ao(naos)) 268 ALLOCATE (almo_scf_env%first_atom_of_domain(ndomains)) 269 ALLOCATE (almo_scf_env%last_atom_of_domain(ndomains)) 270 ALLOCATE (almo_scf_env%nbasis_of_domain(ndomains)) 271 ALLOCATE (almo_scf_env%nocc_of_domain(ndomains, nspins)) !! with smearing, nb of available orbitals for occupation 272 ALLOCATE (almo_scf_env%real_ne_of_domain(ndomains, nspins)) !! with smearing, nb of fully-occupied orbitals 273 ALLOCATE (almo_scf_env%nvirt_full_of_domain(ndomains, nspins)) 274 ALLOCATE (almo_scf_env%nvirt_of_domain(ndomains, nspins)) 275 ALLOCATE (almo_scf_env%nvirt_disc_of_domain(ndomains, nspins)) 276 ALLOCATE (almo_scf_env%mu_of_domain(ndomains, nspins)) 277 ALLOCATE (almo_scf_env%cpu_of_domain(ndomains)) 278 ALLOCATE (almo_scf_env%charge_of_domain(ndomains)) 279 ALLOCATE (almo_scf_env%multiplicity_of_domain(ndomains)) 280 281 ! fill out domain descriptors and group descriptors 282 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular) THEN 283 ! get domain info from molecule_set 284 CALL get_molecule_set_info(molecule_set, & 285 atom_to_mol=almo_scf_env%domain_index_of_atom, & 286 mol_to_first_atom=almo_scf_env%first_atom_of_domain, & 287 mol_to_last_atom=almo_scf_env%last_atom_of_domain, & 288 mol_to_nelectrons=almo_scf_env%nocc_of_domain(1:ndomains, 1), & 289 mol_to_nbasis=almo_scf_env%nbasis_of_domain, & 290 mol_to_charge=almo_scf_env%charge_of_domain, & 291 mol_to_multiplicity=almo_scf_env%multiplicity_of_domain) 292 ! calculate number of alpha and beta occupied orbitals from 293 ! the number of electrons and multiplicity of each molecule 294 ! Na + Nb = Ne 295 ! Na - Nb = Mult - 1 (assume Na > Nb as we do not have more info from get_molecule_set_info) 296 DO idomain = 1, ndomains 297 nelec = almo_scf_env%nocc_of_domain(idomain, 1) 298 multip = almo_scf_env%multiplicity_of_domain(idomain) 299 nelec_a = (nelec + multip - 1)/2 300 nelec_b = nelec - nelec_a 301 !! Initializing an occupation-rescaling trick if smearing is on 302 IF (almo_scf_env%smear) THEN 303 IF (multip .GT. 1) THEN 304 CPWARN("BEWARE: Non singlet state detected, treating it as closed-shell") 305 ENDIF 306 !! Save real number of electrons of each spin, as it is required for Fermi-dirac smearing 307 !! BEWARE : Non singlet states are allowed but treated as closed-shell 308 almo_scf_env%real_ne_of_domain(idomain, :) = REAL(nelec, KIND=dp)/2.0_dp 309 !! Add a number of added_mos equal to the number of atoms in domain 310 !! (since fragments were computed this way with smearing) 311 almo_scf_env%nocc_of_domain(idomain, :) = CEILING(almo_scf_env%real_ne_of_domain(idomain, :)) & 312 + (almo_scf_env%last_atom_of_domain(idomain) & 313 - almo_scf_env%first_atom_of_domain(idomain) + 1) 314 ELSE 315 almo_scf_env%nocc_of_domain(idomain, 1) = nelec_a 316 IF (nelec_a .NE. nelec_b) THEN 317 IF (nspins .EQ. 1) THEN 318 WRITE (*, *) "Domain ", idomain, " out of ", ndomains, ". Electrons = ", nelec 319 CPABORT("odd e- -- use unrestricted methods") 320 ENDIF 321 almo_scf_env%nocc_of_domain(idomain, 2) = nelec_b 322 ! RZK-warning: open-shell procedures have not been tested yet 323 ! Stop the program now 324 CPABORT("Unrestricted ALMO methods are NYI") 325 ENDIF 326 END IF 327 ENDDO 328 DO ispin = 1, nspins 329 ! take care of the full virtual subspace 330 almo_scf_env%nvirt_full_of_domain(:, ispin) = & 331 almo_scf_env%nbasis_of_domain(:) - & 332 almo_scf_env%nocc_of_domain(:, ispin) 333 ! and the truncated virtual subspace 334 SELECT CASE (almo_scf_env%deloc_truncate_virt) 335 CASE (virt_full) 336 almo_scf_env%nvirt_of_domain(:, ispin) = & 337 almo_scf_env%nvirt_full_of_domain(:, ispin) 338 almo_scf_env%nvirt_disc_of_domain(:, ispin) = 0 339 CASE (virt_number) 340 DO idomain = 1, ndomains 341 almo_scf_env%nvirt_of_domain(idomain, ispin) = & 342 MIN(almo_scf_env%deloc_virt_per_domain, & 343 almo_scf_env%nvirt_full_of_domain(idomain, ispin)) 344 almo_scf_env%nvirt_disc_of_domain(idomain, ispin) = & 345 almo_scf_env%nvirt_full_of_domain(idomain, ispin) - & 346 almo_scf_env%nvirt_of_domain(idomain, ispin) 347 ENDDO 348 CASE (virt_occ_size) 349 DO idomain = 1, ndomains 350 almo_scf_env%nvirt_of_domain(idomain, ispin) = & 351 MIN(almo_scf_env%nocc_of_domain(idomain, ispin), & 352 almo_scf_env%nvirt_full_of_domain(idomain, ispin)) 353 almo_scf_env%nvirt_disc_of_domain(idomain, ispin) = & 354 almo_scf_env%nvirt_full_of_domain(idomain, ispin) - & 355 almo_scf_env%nvirt_of_domain(idomain, ispin) 356 ENDDO 357 CASE DEFAULT 358 CPABORT("illegal method for virtual space truncation") 359 END SELECT 360 ENDDO ! spin 361 ELSE ! domains are atomic 362 ! RZK-warning do the same for atomic domains/groups 363 almo_scf_env%domain_index_of_atom(1:natoms) = (/(i, i=1, natoms)/) 364 ENDIF 365 366 ao = 1 367 DO idomain = 1, ndomains 368 DO iao = 1, almo_scf_env%nbasis_of_domain(idomain) 369 almo_scf_env%domain_index_of_ao(ao) = idomain 370 ao = ao + 1 371 ENDDO 372 ENDDO 373 374 almo_scf_env%mu_of_domain(:, :) = almo_scf_env%mu 375 376 ! build domain (i.e. layout) indices for distribution blocks 377 ! ao blocks 378 IF (almo_scf_env%mat_distr_aos == almo_mat_distr_atomic) THEN 379 ALLOCATE (almo_scf_env%domain_index_of_ao_block(natoms)) 380 almo_scf_env%domain_index_of_ao_block(:) = & 381 almo_scf_env%domain_index_of_atom(:) 382 ELSE IF (almo_scf_env%mat_distr_aos == almo_mat_distr_molecular) THEN 383 ALLOCATE (almo_scf_env%domain_index_of_ao_block(nmols)) 384 ! if distr blocks are molecular then domain layout is also molecular 385 almo_scf_env%domain_index_of_ao_block(:) = (/(i, i=1, nmols)/) 386 ENDIF 387 ! mo blocks 388 IF (almo_scf_env%mat_distr_mos == almo_mat_distr_atomic) THEN 389 ALLOCATE (almo_scf_env%domain_index_of_mo_block(natoms)) 390 almo_scf_env%domain_index_of_mo_block(:) = & 391 almo_scf_env%domain_index_of_atom(:) 392 ELSE IF (almo_scf_env%mat_distr_mos == almo_mat_distr_molecular) THEN 393 ALLOCATE (almo_scf_env%domain_index_of_mo_block(nmols)) 394 ! if distr blocks are molecular then domain layout is also molecular 395 almo_scf_env%domain_index_of_mo_block(:) = (/(i, i=1, nmols)/) 396 ENDIF 397 398 ! set all flags 399 !almo_scf_env%need_previous_ks=.FALSE. 400 !IF (almo_scf_env%deloc_method==almo_deloc_harris) THEN 401 almo_scf_env%need_previous_ks = .TRUE. 402 !ENDIF 403 404 !almo_scf_env%need_virtuals=.FALSE. 405 !almo_scf_env%need_orbital_energies=.FALSE. 406 !IF (almo_scf_env%almo_update_algorithm==almo_scf_diag) THEN 407 almo_scf_env%need_virtuals = .TRUE. 408 almo_scf_env%need_orbital_energies = .TRUE. 409 !ENDIF 410 411 almo_scf_env%calc_forces = calc_forces 412 IF (calc_forces) THEN 413 CALL cite_reference(Scheiber2018) 414 IF (almo_scf_env%deloc_method .EQ. almo_deloc_x .OR. & 415 almo_scf_env%deloc_method .EQ. almo_deloc_xalmo_x .OR. & 416 almo_scf_env%deloc_method .EQ. almo_deloc_xalmo_1diag) THEN 417 CPABORT("Forces for perturbative methods are NYI. Change DELOCALIZE_METHOD") 418 ENDIF 419 ! switch to ASPC after a certain number of exact steps is done 420 IF (almo_scf_env%almo_history%istore .GT. (almo_scf_env%almo_history%nstore + 1)) THEN 421 IF (almo_scf_env%opt_block_diag_pcg%eps_error_early .GT. 0.0_dp) THEN 422 almo_scf_env%opt_block_diag_pcg%eps_error = almo_scf_env%opt_block_diag_pcg%eps_error_early 423 almo_scf_env%opt_block_diag_pcg%early_stopping_on = .TRUE. 424 IF (unit_nr > 0) WRITE (*, *) "ALMO_OPTIMIZER_PCG: EPS_ERROR_EARLY is on" 425 ENDIF 426 IF (almo_scf_env%opt_block_diag_diis%eps_error_early .GT. 0.0_dp) THEN 427 almo_scf_env%opt_block_diag_diis%eps_error = almo_scf_env%opt_block_diag_diis%eps_error_early 428 almo_scf_env%opt_block_diag_diis%early_stopping_on = .TRUE. 429 IF (unit_nr > 0) WRITE (*, *) "ALMO_OPTIMIZER_DIIS: EPS_ERROR_EARLY is on" 430 ENDIF 431 IF (almo_scf_env%opt_block_diag_pcg%max_iter_early .GT. 0) THEN 432 almo_scf_env%opt_block_diag_pcg%max_iter = almo_scf_env%opt_block_diag_pcg%max_iter_early 433 almo_scf_env%opt_block_diag_pcg%early_stopping_on = .TRUE. 434 IF (unit_nr > 0) WRITE (*, *) "ALMO_OPTIMIZER_PCG: MAX_ITER_EARLY is on" 435 ENDIF 436 IF (almo_scf_env%opt_block_diag_diis%max_iter_early .GT. 0) THEN 437 almo_scf_env%opt_block_diag_diis%max_iter = almo_scf_env%opt_block_diag_diis%max_iter_early 438 almo_scf_env%opt_block_diag_diis%early_stopping_on = .TRUE. 439 IF (unit_nr > 0) WRITE (*, *) "ALMO_OPTIMIZER_DIIS: MAX_ITER_EARLY is on" 440 ENDIF 441 ELSE 442 almo_scf_env%opt_block_diag_diis%early_stopping_on = .FALSE. 443 almo_scf_env%opt_block_diag_pcg%early_stopping_on = .FALSE. 444 ENDIF 445 IF (almo_scf_env%xalmo_history%istore .GT. (almo_scf_env%xalmo_history%nstore + 1)) THEN 446 IF (almo_scf_env%opt_xalmo_pcg%eps_error_early .GT. 0.0_dp) THEN 447 almo_scf_env%opt_xalmo_pcg%eps_error = almo_scf_env%opt_xalmo_pcg%eps_error_early 448 almo_scf_env%opt_xalmo_pcg%early_stopping_on = .TRUE. 449 IF (unit_nr > 0) WRITE (*, *) "XALMO_OPTIMIZER_PCG: EPS_ERROR_EARLY is on" 450 ENDIF 451 IF (almo_scf_env%opt_xalmo_pcg%max_iter_early .GT. 0.0_dp) THEN 452 almo_scf_env%opt_xalmo_pcg%max_iter = almo_scf_env%opt_xalmo_pcg%max_iter_early 453 almo_scf_env%opt_xalmo_pcg%early_stopping_on = .TRUE. 454 IF (unit_nr > 0) WRITE (*, *) "XALMO_OPTIMIZER_PCG: MAX_ITER_EARLY is on" 455 ENDIF 456 ELSE 457 almo_scf_env%opt_xalmo_pcg%early_stopping_on = .FALSE. 458 ENDIF 459 ENDIF 460 461 ! create all matrices 462 CALL almo_scf_env_create_matrices(almo_scf_env, matrix_s(1)%matrix) 463 464 ! set up matrix S and all required functions of S 465 almo_scf_env%s_inv_done = .FALSE. 466 almo_scf_env%s_sqrt_done = .FALSE. 467 CALL almo_scf_init_ao_overlap(matrix_s(1)%matrix, almo_scf_env) 468 469 ! create the quencher (imposes sparsity template) 470 CALL almo_scf_construct_quencher(qs_env, almo_scf_env) 471 CALL distribute_domains(almo_scf_env) 472 473 ! FINISH setting job parameters here, print out job info 474 CALL almo_scf_print_job_info(almo_scf_env, unit_nr) 475 476 ! allocate and init the domain preconditioner 477 ALLOCATE (almo_scf_env%domain_preconditioner(ndomains, nspins)) 478 CALL init_submatrices(almo_scf_env%domain_preconditioner) 479 480 ! allocate and init projected KS for domains 481 ALLOCATE (almo_scf_env%domain_ks_xx(ndomains, nspins)) 482 CALL init_submatrices(almo_scf_env%domain_ks_xx) 483 484 ! init ao-overlap subblocks 485 ALLOCATE (almo_scf_env%domain_s_inv(ndomains, nspins)) 486 CALL init_submatrices(almo_scf_env%domain_s_inv) 487 ALLOCATE (almo_scf_env%domain_s_sqrt_inv(ndomains, nspins)) 488 CALL init_submatrices(almo_scf_env%domain_s_sqrt_inv) 489 ALLOCATE (almo_scf_env%domain_s_sqrt(ndomains, nspins)) 490 CALL init_submatrices(almo_scf_env%domain_s_sqrt) 491 ALLOCATE (almo_scf_env%domain_t(ndomains, nspins)) 492 CALL init_submatrices(almo_scf_env%domain_t) 493 ALLOCATE (almo_scf_env%domain_err(ndomains, nspins)) 494 CALL init_submatrices(almo_scf_env%domain_err) 495 ALLOCATE (almo_scf_env%domain_r_down_up(ndomains, nspins)) 496 CALL init_submatrices(almo_scf_env%domain_r_down_up) 497 498 ! initialization of the KS matrix 499 CALL init_almo_ks_matrix_via_qs(qs_env, & 500 almo_scf_env%matrix_ks, & 501 almo_scf_env%mat_distr_aos, & 502 almo_scf_env%eps_filter) 503 CALL construct_qs_mos(qs_env, almo_scf_env) 504 505 CALL timestop(handle) 506 507 END SUBROUTINE almo_scf_init 508 509! ************************************************************************************************** 510!> \brief create the scf initial guess for ALMOs 511!> \param qs_env ... 512!> \param almo_scf_env ... 513!> \par History 514!> 2016.11 created [Rustam Z Khaliullin] 515!> 2018.09 smearing support [Ruben Staub] 516!> \author Rustam Z Khaliullin 517! ************************************************************************************************** 518 SUBROUTINE almo_scf_initial_guess(qs_env, almo_scf_env) 519 TYPE(qs_environment_type), POINTER :: qs_env 520 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 521 522 CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_initial_guess', & 523 routineP = moduleN//':'//routineN 524 525 CHARACTER(LEN=default_path_length) :: file_name, project_name 526 INTEGER :: handle, iaspc, ispin, istore, naspc, & 527 nspins, unit_nr 528 INTEGER, DIMENSION(2) :: nelectron_spin 529 LOGICAL :: aspc_guess, has_unit_metric 530 REAL(KIND=dp) :: alpha, cs_pos, energy, kTS_sum 531 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set 532 TYPE(cp_logger_type), POINTER :: logger 533 TYPE(cp_para_env_type), POINTER :: para_env 534 TYPE(dbcsr_distribution_type) :: dist 535 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s, rho_ao 536 TYPE(dft_control_type), POINTER :: dft_control 537 TYPE(molecular_scf_guess_env_type), POINTER :: mscfg_env 538 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set 539 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set 540 TYPE(qs_rho_type), POINTER :: rho 541 542 CALL timeset(routineN, handle) 543 544 NULLIFY (rho, rho_ao) 545 546 ! get a useful output_unit 547 logger => cp_get_default_logger() 548 IF (logger%para_env%ionode) THEN 549 unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.) 550 ELSE 551 unit_nr = -1 552 ENDIF 553 554 ! get basic quantities from the qs_env 555 CALL get_qs_env(qs_env, & 556 dft_control=dft_control, & 557 matrix_s=matrix_s, & 558 atomic_kind_set=atomic_kind_set, & 559 qs_kind_set=qs_kind_set, & 560 particle_set=particle_set, & 561 has_unit_metric=has_unit_metric, & 562 para_env=para_env, & 563 nelectron_spin=nelectron_spin, & 564 mscfg_env=mscfg_env, & 565 rho=rho) 566 567 CALL qs_rho_get(rho, rho_ao=rho_ao) 568 CPASSERT(ASSOCIATED(mscfg_env)) 569 570 ! initial guess on the first simulation step is determined by almo_scf_env%almo_scf_guess 571 ! the subsequent simulation steps are determined by extrapolation_order 572 ! if extrapolation order is zero then again almo_scf_env%almo_scf_guess is used 573 ! ... the number of stored history points will remain zero if extrapolation order is zero 574 IF (almo_scf_env%almo_history%istore == 0) THEN 575 aspc_guess = .FALSE. 576 ELSE 577 aspc_guess = .TRUE. 578 ENDIF 579 580 nspins = almo_scf_env%nspins 581 582 ! create an initial guess 583 IF (.NOT. aspc_guess) THEN 584 585 SELECT CASE (almo_scf_env%almo_scf_guess) 586 CASE (molecular_guess) 587 588 DO ispin = 1, nspins 589 590 ! the calculations on "isolated" molecules has already been done 591 ! all we need to do is convert the MOs of molecules into 592 ! the ALMO matrix taking into account different distributions 593 CALL get_matrix_from_submatrices(mscfg_env, & 594 almo_scf_env%matrix_t_blk(ispin), ispin) 595 CALL dbcsr_filter(almo_scf_env%matrix_t_blk(ispin), & 596 almo_scf_env%eps_filter) 597 598 ENDDO 599 600 CASE (atomic_guess) 601 602 IF (dft_control%qs_control%dftb .OR. dft_control%qs_control%semi_empirical .OR. & 603 dft_control%qs_control%xtb) THEN 604 CALL calculate_mopac_dm(rho_ao, & 605 matrix_s(1)%matrix, has_unit_metric, & 606 dft_control, particle_set, atomic_kind_set, qs_kind_set, & 607 nspins, nelectron_spin, & 608 para_env) 609 ELSE 610 CALL calculate_atomic_block_dm(rho_ao, & 611 matrix_s(1)%matrix, particle_set, atomic_kind_set, qs_kind_set, & 612 nspins, nelectron_spin, unit_nr, para_env) 613 ENDIF 614 615 DO ispin = 1, nspins 616 ! copy the atomic-block dm into matrix_p_blk 617 CALL matrix_qs_to_almo(rho_ao(ispin)%matrix, & 618 almo_scf_env%matrix_p_blk(ispin), almo_scf_env%mat_distr_aos, & 619 .FALSE.) 620 CALL dbcsr_filter(almo_scf_env%matrix_p_blk(ispin), & 621 almo_scf_env%eps_filter) 622 ENDDO ! ispin 623 624 ! obtain orbitals from the density matrix 625 ! (the current version of ALMO SCF needs orbitals) 626 CALL almo_scf_p_blk_to_t_blk(almo_scf_env, ionic=.FALSE.) 627 628 CASE (restart_guess) 629 630 project_name = logger%iter_info%project_name 631 632 DO ispin = 1, nspins 633 WRITE (file_name, '(A,I0,A)') TRIM(project_name)//"_ALMO_SPIN_", ispin, "_RESTART.mo" 634 CALL dbcsr_get_info(almo_scf_env%matrix_t_blk(ispin), distribution=dist) 635 CALL dbcsr_binary_read(file_name, distribution=dist, matrix_new=almo_scf_env%matrix_t_blk(ispin)) 636 cs_pos = dbcsr_checksum(almo_scf_env%matrix_t_blk(ispin), pos=.TRUE.) 637 IF (unit_nr > 0) THEN 638 WRITE (unit_nr, '(T2,A,E20.8)') "Read restart ALMO "//TRIM(file_name)//" with checksum: ", cs_pos 639 ENDIF 640 ENDDO 641 END SELECT 642 643 ELSE !aspc_guess 644 645 CALL cite_reference(Kolafa2004) 646 647 naspc = MIN(almo_scf_env%almo_history%istore, almo_scf_env%almo_history%nstore) 648 IF (unit_nr > 0) THEN 649 WRITE (unit_nr, FMT="(/,T2,A,/,/,T3,A,I0)") & 650 "Parameters for the always stable predictor-corrector (ASPC) method:", & 651 "ASPC order: ", naspc 652 END IF 653 654 DO ispin = 1, nspins 655 656 ! extrapolation 657 DO iaspc = 1, naspc 658 istore = MOD(almo_scf_env%almo_history%istore - iaspc, almo_scf_env%almo_history%nstore) + 1 659 alpha = (-1.0_dp)**(iaspc + 1)*REAL(iaspc, KIND=dp)* & 660 binomial(2*naspc, naspc - iaspc)/binomial(2*naspc - 2, naspc - 1) 661 IF (unit_nr > 0) THEN 662 WRITE (unit_nr, FMT="(T3,A2,I0,A4,F10.6)") & 663 "B(", iaspc, ") = ", alpha 664 END IF 665 IF (iaspc == 1) THEN 666 CALL dbcsr_copy(almo_scf_env%matrix_t_blk(ispin), & 667 almo_scf_env%almo_history%matrix_t(ispin), & 668 keep_sparsity=.TRUE.) 669 CALL dbcsr_scale(almo_scf_env%matrix_t_blk(ispin), alpha) 670 ELSE 671 CALL dbcsr_multiply("N", "N", alpha, & 672 almo_scf_env%almo_history%matrix_p_up_down(ispin, istore), & 673 almo_scf_env%almo_history%matrix_t(ispin), & 674 1.0_dp, almo_scf_env%matrix_t_blk(ispin), & 675 retain_sparsity=.TRUE.) 676 ENDIF 677 ENDDO !iaspc 678 679 ENDDO !ispin 680 681 ENDIF !aspc_guess? 682 683 DO ispin = 1, nspins 684 685 CALL orthogonalize_mos(ket=almo_scf_env%matrix_t_blk(ispin), & 686 overlap=almo_scf_env%matrix_sigma_blk(ispin), & 687 metric=almo_scf_env%matrix_s_blk(1), & 688 retain_locality=.TRUE., & 689 only_normalize=.FALSE., & 690 nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), & 691 eps_filter=almo_scf_env%eps_filter, & 692 order_lanczos=almo_scf_env%order_lanczos, & 693 eps_lanczos=almo_scf_env%eps_lanczos, & 694 max_iter_lanczos=almo_scf_env%max_iter_lanczos) 695 696 !! Application of an occupation-rescaling trick for smearing, if requested 697 IF (almo_scf_env%smear) THEN 698 CALL almo_scf_t_rescaling(matrix_t=almo_scf_env%matrix_t_blk(ispin), & 699 mo_energies=almo_scf_env%mo_energies(:, ispin), & 700 mu_of_domain=almo_scf_env%mu_of_domain(:, ispin), & 701 real_ne_of_domain=almo_scf_env%real_ne_of_domain(:, ispin), & 702 spin_kTS=almo_scf_env%kTS(ispin), & 703 smear_e_temp=almo_scf_env%smear_e_temp, & 704 ndomains=almo_scf_env%ndomains, & 705 nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin)) 706 END IF 707 708 CALL almo_scf_t_to_proj(t=almo_scf_env%matrix_t_blk(ispin), & 709 p=almo_scf_env%matrix_p(ispin), & 710 eps_filter=almo_scf_env%eps_filter, & 711 orthog_orbs=.FALSE., & 712 nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), & 713 s=almo_scf_env%matrix_s(1), & 714 sigma=almo_scf_env%matrix_sigma(ispin), & 715 sigma_inv=almo_scf_env%matrix_sigma_inv(ispin), & 716 use_guess=.FALSE., & 717 smear=almo_scf_env%smear, & 718 algorithm=almo_scf_env%sigma_inv_algorithm, & 719 eps_lanczos=almo_scf_env%eps_lanczos, & 720 max_iter_lanczos=almo_scf_env%max_iter_lanczos, & 721 inv_eps_factor=almo_scf_env%matrix_iter_eps_error_factor, & 722 para_env=almo_scf_env%para_env, & 723 blacs_env=almo_scf_env%blacs_env) 724 725 ENDDO 726 727 ! compute dm from the projector(s) 728 IF (nspins == 1) THEN 729 CALL dbcsr_scale(almo_scf_env%matrix_p(1), 2.0_dp) 730 !! Rescaling electronic entropy contribution by spin_factor 731 IF (almo_scf_env%smear) THEN 732 almo_scf_env%kTS(1) = almo_scf_env%kTS(1)*2.0_dp 733 END IF 734 ENDIF 735 736 IF (almo_scf_env%smear) THEN 737 kTS_sum = SUM(almo_scf_env%kTS) 738 ELSE 739 kTS_sum = 0.0_dp 740 ENDIF 741 742 CALL almo_dm_to_almo_ks(qs_env, & 743 almo_scf_env%matrix_p, & 744 almo_scf_env%matrix_ks, & 745 energy, & 746 almo_scf_env%eps_filter, & 747 almo_scf_env%mat_distr_aos, & 748 smear=almo_scf_env%smear, & 749 kTS_sum=kTS_sum) 750 751 IF (unit_nr > 0) THEN 752 IF (almo_scf_env%almo_scf_guess .EQ. molecular_guess) THEN 753 WRITE (unit_nr, '(T2,A38,F40.10)') "Single-molecule energy:", & 754 SUM(mscfg_env%energy_of_frag) 755 ENDIF 756 WRITE (unit_nr, '(T2,A38,F40.10)') "Energy of the initial guess:", energy 757 WRITE (unit_nr, '()') 758 ENDIF 759 760 CALL timestop(handle) 761 762 END SUBROUTINE almo_scf_initial_guess 763 764! ************************************************************************************************** 765!> \brief store a history of matrices for later use in almo_scf_initial_guess 766!> \param almo_scf_env ... 767!> \par History 768!> 2016.11 created [Rustam Z Khaliullin] 769!> \author Rustam Khaliullin 770! ************************************************************************************************** 771 SUBROUTINE almo_scf_store_extrapolation_data(almo_scf_env) 772 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 773 774 CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_store_extrapolation_data', & 775 routineP = moduleN//':'//routineN 776 777 INTEGER :: handle, ispin, istore, unit_nr 778 LOGICAL :: delocalization_uses_extrapolation 779 TYPE(cp_logger_type), POINTER :: logger 780 TYPE(dbcsr_type) :: matrix_no_tmp1, matrix_no_tmp2, & 781 matrix_no_tmp3, matrix_no_tmp4 782 783 CALL timeset(routineN, handle) 784 785 ! get a useful output_unit 786 logger => cp_get_default_logger() 787 IF (logger%para_env%ionode) THEN 788 unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.) 789 ELSE 790 unit_nr = -1 791 ENDIF 792 793 IF (almo_scf_env%almo_history%nstore > 0) THEN 794 795 almo_scf_env%almo_history%istore = almo_scf_env%almo_history%istore + 1 796 797 DO ispin = 1, SIZE(almo_scf_env%matrix_t_blk) 798 799 istore = MOD(almo_scf_env%almo_history%istore - 1, almo_scf_env%almo_history%nstore) + 1 800 801 IF (almo_scf_env%almo_history%istore == 1) THEN 802 CALL dbcsr_create(almo_scf_env%almo_history%matrix_t(ispin), & 803 template=almo_scf_env%matrix_t_blk(ispin), & 804 matrix_type=dbcsr_type_no_symmetry) 805 ENDIF 806 CALL dbcsr_copy(almo_scf_env%almo_history%matrix_t(ispin), & 807 almo_scf_env%matrix_t_blk(ispin)) 808 809 IF (almo_scf_env%almo_history%istore <= almo_scf_env%almo_history%nstore) THEN 810 CALL dbcsr_create(almo_scf_env%almo_history%matrix_p_up_down(ispin, istore), & 811 template=almo_scf_env%matrix_s(1), & 812 matrix_type=dbcsr_type_no_symmetry) 813 ENDIF 814 815 CALL dbcsr_create(matrix_no_tmp1, template=almo_scf_env%matrix_t_blk(ispin), & 816 matrix_type=dbcsr_type_no_symmetry) 817 CALL dbcsr_create(matrix_no_tmp2, template=almo_scf_env%matrix_t_blk(ispin), & 818 matrix_type=dbcsr_type_no_symmetry) 819 820 ! compute contra-covariant density matrix 821 CALL dbcsr_multiply("N", "N", 1.0_dp, almo_scf_env%matrix_s(1), & 822 almo_scf_env%matrix_t_blk(ispin), & 823 0.0_dp, matrix_no_tmp1, & 824 filter_eps=almo_scf_env%eps_filter) 825 CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_no_tmp1, & 826 almo_scf_env%matrix_sigma_inv_0deloc(ispin), & 827 0.0_dp, matrix_no_tmp2, & 828 filter_eps=almo_scf_env%eps_filter) 829 CALL dbcsr_multiply("N", "T", 1.0_dp, & 830 almo_scf_env%matrix_t_blk(ispin), & 831 matrix_no_tmp2, & 832 0.0_dp, almo_scf_env%almo_history%matrix_p_up_down(ispin, istore), & 833 filter_eps=almo_scf_env%eps_filter) 834 835 CALL dbcsr_release(matrix_no_tmp1) 836 CALL dbcsr_release(matrix_no_tmp2) 837 838 ENDDO 839 840 ENDIF 841 842 ! exrapolate xalmos? 843 delocalization_uses_extrapolation = & 844 .NOT. ((almo_scf_env%deloc_method .EQ. almo_deloc_none) .OR. & 845 (almo_scf_env%deloc_method .EQ. almo_deloc_xalmo_1diag)) 846 IF (almo_scf_env%xalmo_history%nstore > 0 .AND. & 847 delocalization_uses_extrapolation) THEN 848 849 almo_scf_env%xalmo_history%istore = almo_scf_env%xalmo_history%istore + 1 850 851 DO ispin = 1, SIZE(almo_scf_env%matrix_t) 852 853 istore = MOD(almo_scf_env%xalmo_history%istore - 1, almo_scf_env%xalmo_history%nstore) + 1 854 855 IF (almo_scf_env%xalmo_history%istore == 1) THEN 856 CALL dbcsr_create(almo_scf_env%xalmo_history%matrix_t(ispin), & 857 template=almo_scf_env%matrix_t(ispin), & 858 matrix_type=dbcsr_type_no_symmetry) 859 ENDIF 860 CALL dbcsr_copy(almo_scf_env%xalmo_history%matrix_t(ispin), & 861 almo_scf_env%matrix_t(ispin)) 862 863 IF (almo_scf_env%xalmo_history%istore <= almo_scf_env%xalmo_history%nstore) THEN 864 !CALL dbcsr_init(almo_scf_env%xalmo_history%matrix_x(ispin, istore)) 865 !CALL dbcsr_create(almo_scf_env%xalmo_history%matrix_x(ispin, istore), & 866 ! template=almo_scf_env%matrix_t(ispin), & 867 ! matrix_type=dbcsr_type_no_symmetry) 868 CALL dbcsr_create(almo_scf_env%xalmo_history%matrix_p_up_down(ispin, istore), & 869 template=almo_scf_env%matrix_s(1), & 870 matrix_type=dbcsr_type_no_symmetry) 871 ENDIF 872 873 CALL dbcsr_create(matrix_no_tmp3, template=almo_scf_env%matrix_t(ispin), & 874 matrix_type=dbcsr_type_no_symmetry) 875 CALL dbcsr_create(matrix_no_tmp4, template=almo_scf_env%matrix_t(ispin), & 876 matrix_type=dbcsr_type_no_symmetry) 877 878 ! compute contra-covariant density matrix 879 CALL dbcsr_multiply("N", "N", 1.0_dp, almo_scf_env%matrix_s(1), & 880 almo_scf_env%matrix_t(ispin), & 881 0.0_dp, matrix_no_tmp3, & 882 filter_eps=almo_scf_env%eps_filter) 883 CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_no_tmp3, & 884 almo_scf_env%matrix_sigma_inv(ispin), & 885 0.0_dp, matrix_no_tmp4, & 886 filter_eps=almo_scf_env%eps_filter) 887 CALL dbcsr_multiply("N", "T", 1.0_dp, & 888 almo_scf_env%matrix_t(ispin), & 889 matrix_no_tmp4, & 890 0.0_dp, almo_scf_env%xalmo_history%matrix_p_up_down(ispin, istore), & 891 filter_eps=almo_scf_env%eps_filter) 892 893 ! store the difference between t and t0 894 !CALL dbcsr_copy(almo_scf_env%xalmo_history%matrix_x(ispin, istore),& 895 ! almo_scf_env%matrix_t(ispin)) 896 !CALL dbcsr_add(almo_scf_env%xalmo_history%matrix_x(ispin, istore),& 897 ! almo_scf_env%matrix_t_blk(ispin),1.0_dp,-1.0_dp) 898 899 CALL dbcsr_release(matrix_no_tmp3) 900 CALL dbcsr_release(matrix_no_tmp4) 901 902 ENDDO 903 904 ENDIF 905 906 CALL timestop(handle) 907 908 END SUBROUTINE almo_scf_store_extrapolation_data 909 910! ************************************************************************************************** 911!> \brief Prints out a short summary about the ALMO SCF job 912!> \param almo_scf_env ... 913!> \param unit_nr ... 914!> \par History 915!> 2011.10 created [Rustam Z Khaliullin] 916!> \author Rustam Z Khaliullin 917! ************************************************************************************************** 918 SUBROUTINE almo_scf_print_job_info(almo_scf_env, unit_nr) 919 920 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 921 INTEGER, INTENT(IN) :: unit_nr 922 923 CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_print_job_info', & 924 routineP = moduleN//':'//routineN 925 926 CHARACTER(len=13) :: neig_string 927 CHARACTER(len=33) :: deloc_method_string 928 INTEGER :: handle, idomain, index1_prev, sum_temp 929 INTEGER, ALLOCATABLE, DIMENSION(:) :: nneighbors 930 931 CALL timeset(routineN, handle) 932 933 IF (unit_nr > 0) THEN 934 WRITE (unit_nr, '()') 935 WRITE (unit_nr, '(T2,A,A,A)') REPEAT("-", 32), " ALMO SETTINGS ", REPEAT("-", 32) 936 937 WRITE (unit_nr, '(T2,A,T48,E33.3)') "eps_filter:", almo_scf_env%eps_filter 938 939 IF (almo_scf_env%almo_update_algorithm .EQ. almo_scf_skip) THEN 940 WRITE (unit_nr, '(T2,A)') "skip optimization of block-diagonal ALMOs" 941 ELSE 942 WRITE (unit_nr, '(T2,A)') "optimization of block-diagonal ALMOs:" 943 SELECT CASE (almo_scf_env%almo_update_algorithm) 944 CASE (almo_scf_diag) 945 ! the DIIS algorith is the only choice for the diagonlaization-based algorithm 946 CALL print_optimizer_options(almo_scf_env%opt_block_diag_diis, unit_nr) 947 CASE (almo_scf_pcg) 948 ! print out PCG options 949 CALL print_optimizer_options(almo_scf_env%opt_block_diag_pcg, unit_nr) 950 END SELECT 951 ENDIF 952 953 SELECT CASE (almo_scf_env%deloc_method) 954 CASE (almo_deloc_none) 955 deloc_method_string = "NONE" 956 CASE (almo_deloc_x) 957 deloc_method_string = "FULL_X" 958 CASE (almo_deloc_scf) 959 deloc_method_string = "FULL_SCF" 960 CASE (almo_deloc_x_then_scf) 961 deloc_method_string = "FULL_X_THEN_SCF" 962 CASE (almo_deloc_xalmo_1diag) 963 deloc_method_string = "XALMO_1DIAG" 964 CASE (almo_deloc_xalmo_x) 965 deloc_method_string = "XALMO_X" 966 CASE (almo_deloc_xalmo_scf) 967 deloc_method_string = "XALMO_SCF" 968 END SELECT 969 WRITE (unit_nr, '(T2,A,T48,A33)') "delocalization:", TRIM(deloc_method_string) 970 971 IF (almo_scf_env%deloc_method .NE. almo_deloc_none) THEN 972 973 SELECT CASE (almo_scf_env%deloc_method) 974 CASE (almo_deloc_x, almo_deloc_scf, almo_deloc_x_then_scf) 975 WRITE (unit_nr, '(T2,A,T48,A33)') "delocalization cutoff radius:", & 976 "infinite" 977 deloc_method_string = "FULL_X_THEN_SCF" 978 CASE (almo_deloc_xalmo_1diag, almo_deloc_xalmo_x, almo_deloc_xalmo_scf) 979 WRITE (unit_nr, '(T2,A,T48,F33.5)') "XALMO cutoff radius:", & 980 almo_scf_env%quencher_r0_factor 981 END SELECT 982 983 IF (almo_scf_env%deloc_method .EQ. almo_deloc_xalmo_1diag) THEN 984 ! print nothing because no actual optimization is done 985 ELSE 986 WRITE (unit_nr, '(T2,A)') "optimization of extended orbitals:" 987 SELECT CASE (almo_scf_env%xalmo_update_algorithm) 988 CASE (almo_scf_diag) 989 CALL print_optimizer_options(almo_scf_env%opt_xalmo_diis, unit_nr) 990 CASE (almo_scf_trustr) 991 CALL print_optimizer_options(almo_scf_env%opt_xalmo_trustr, unit_nr) 992 CASE (almo_scf_pcg) 993 CALL print_optimizer_options(almo_scf_env%opt_xalmo_pcg, unit_nr) 994 END SELECT 995 ENDIF 996 997 ENDIF 998 999 !SELECT CASE(almo_scf_env%domain_layout_mos) 1000 !CASE(almo_domain_layout_orbital) 1001 ! WRITE(unit_nr,'(T2,A,T48,A33)') "Delocalization constraints","ORBITAL" 1002 !CASE(almo_domain_layout_atomic) 1003 ! WRITE(unit_nr,'(T2,A,T48,A33)') "Delocalization constraints","ATOMIC" 1004 !CASE(almo_domain_layout_molecular) 1005 ! WRITE(unit_nr,'(T2,A,T48,A33)') "Delocalization constraints","MOLECULAR" 1006 !END SELECT 1007 1008 !SELECT CASE(almo_scf_env%domain_layout_aos) 1009 !CASE(almo_domain_layout_atomic) 1010 ! WRITE(unit_nr,'(T2,A,T48,A33)') "Basis function domains","ATOMIC" 1011 !CASE(almo_domain_layout_molecular) 1012 ! WRITE(unit_nr,'(T2,A,T48,A33)') "Basis function domains","MOLECULAR" 1013 !END SELECT 1014 1015 !SELECT CASE(almo_scf_env%mat_distr_aos) 1016 !CASE(almo_mat_distr_atomic) 1017 ! WRITE(unit_nr,'(T2,A,T48,A33)') "Parallel distribution for AOs","ATOMIC" 1018 !CASE(almo_mat_distr_molecular) 1019 ! WRITE(unit_nr,'(T2,A,T48,A33)') "Parallel distribution for AOs","MOLECULAR" 1020 !END SELECT 1021 1022 !SELECT CASE(almo_scf_env%mat_distr_mos) 1023 !CASE(almo_mat_distr_atomic) 1024 ! WRITE(unit_nr,'(T2,A,T48,A33)') "Parallel distribution for MOs","ATOMIC" 1025 !CASE(almo_mat_distr_molecular) 1026 ! WRITE(unit_nr,'(T2,A,T48,A33)') "Parallel distribution for MOs","MOLECULAR" 1027 !END SELECT 1028 1029 ! print fragment's statistics 1030 WRITE (unit_nr, '(T2,A)') REPEAT("-", 79) 1031 WRITE (unit_nr, '(T2,A,T48,I33)') "Total fragments:", & 1032 almo_scf_env%ndomains 1033 1034 sum_temp = SUM(almo_scf_env%nbasis_of_domain(:)) 1035 WRITE (unit_nr, '(T2,A,T53,I5,F9.2,I5,I9)') & 1036 "Basis set size per fragment (min, av, max, total):", & 1037 MINVAL(almo_scf_env%nbasis_of_domain(:)), & 1038 (1.0_dp*sum_temp)/almo_scf_env%ndomains, & 1039 MAXVAL(almo_scf_env%nbasis_of_domain(:)), & 1040 sum_temp 1041 !WRITE (unit_nr, '(T2,I13,F13.3,I13,I13)') & 1042 ! MINVAL(almo_scf_env%nbasis_of_domain(:)), & 1043 ! (1.0_dp*sum_temp) / almo_scf_env%ndomains, & 1044 ! MAXVAL(almo_scf_env%nbasis_of_domain(:)), & 1045 ! sum_temp 1046 1047 sum_temp = SUM(almo_scf_env%nocc_of_domain(:, :)) 1048 WRITE (unit_nr, '(T2,A,T53,I5,F9.2,I5,I9)') & 1049 "Occupied MOs per fragment (min, av, max, total):", & 1050 MINVAL(SUM(almo_scf_env%nocc_of_domain, DIM=2)), & 1051 (1.0_dp*sum_temp)/almo_scf_env%ndomains, & 1052 MAXVAL(SUM(almo_scf_env%nocc_of_domain, DIM=2)), & 1053 sum_temp 1054 !WRITE (unit_nr, '(T2,I13,F13.3,I13,I13)') & 1055 ! MINVAL( SUM(almo_scf_env%nocc_of_domain, DIM=2) ), & 1056 ! (1.0_dp*sum_temp) / almo_scf_env%ndomains, & 1057 ! MAXVAL( SUM(almo_scf_env%nocc_of_domain, DIM=2) ), & 1058 ! sum_temp 1059 1060 sum_temp = SUM(almo_scf_env%nvirt_of_domain(:, :)) 1061 WRITE (unit_nr, '(T2,A,T53,I5,F9.2,I5,I9)') & 1062 "Virtual MOs per fragment (min, av, max, total):", & 1063 MINVAL(SUM(almo_scf_env%nvirt_of_domain, DIM=2)), & 1064 (1.0_dp*sum_temp)/almo_scf_env%ndomains, & 1065 MAXVAL(SUM(almo_scf_env%nvirt_of_domain, DIM=2)), & 1066 sum_temp 1067 !WRITE (unit_nr, '(T2,I13,F13.3,I13,I13)') & 1068 ! MINVAL( SUM(almo_scf_env%nvirt_of_domain, DIM=2) ), & 1069 ! (1.0_dp*sum_temp) / almo_scf_env%ndomains, & 1070 ! MAXVAL( SUM(almo_scf_env%nvirt_of_domain, DIM=2) ), & 1071 ! sum_temp 1072 1073 sum_temp = SUM(almo_scf_env%charge_of_domain(:)) 1074 WRITE (unit_nr, '(T2,A,T53,I5,F9.2,I5,I9)') & 1075 "Charges per fragment (min, av, max, total):", & 1076 MINVAL(almo_scf_env%charge_of_domain(:)), & 1077 (1.0_dp*sum_temp)/almo_scf_env%ndomains, & 1078 MAXVAL(almo_scf_env%charge_of_domain(:)), & 1079 sum_temp 1080 !WRITE (unit_nr, '(T2,I13,F13.3,I13,I13)') & 1081 ! MINVAL(almo_scf_env%charge_of_domain(:)), & 1082 ! (1.0_dp*sum_temp) / almo_scf_env%ndomains, & 1083 ! MAXVAL(almo_scf_env%charge_of_domain(:)), & 1084 ! sum_temp 1085 1086 ! compute the number of neighbors of each fragment 1087 ALLOCATE (nneighbors(almo_scf_env%ndomains)) 1088 1089 DO idomain = 1, almo_scf_env%ndomains 1090 1091 IF (idomain .EQ. 1) THEN 1092 index1_prev = 1 1093 ELSE 1094 index1_prev = almo_scf_env%domain_map(1)%index1(idomain - 1) 1095 ENDIF 1096 1097 SELECT CASE (almo_scf_env%deloc_method) 1098 CASE (almo_deloc_none) 1099 nneighbors(idomain) = 0 1100 CASE (almo_deloc_x, almo_deloc_scf, almo_deloc_x_then_scf) 1101 nneighbors(idomain) = almo_scf_env%ndomains - 1 ! minus self 1102 CASE (almo_deloc_xalmo_1diag, almo_deloc_xalmo_x, almo_deloc_xalmo_scf) 1103 nneighbors(idomain) = almo_scf_env%domain_map(1)%index1(idomain) - index1_prev - 1 ! minus self 1104 CASE DEFAULT 1105 nneighbors(idomain) = -1 1106 END SELECT 1107 1108 ENDDO ! cycle over domains 1109 1110 sum_temp = SUM(nneighbors(:)) 1111 WRITE (unit_nr, '(T2,A,T53,I5,F9.2,I5,I9)') & 1112 "Deloc. neighbors of fragment (min, av, max, total):", & 1113 MINVAL(nneighbors(:)), & 1114 (1.0_dp*sum_temp)/almo_scf_env%ndomains, & 1115 MAXVAL(nneighbors(:)), & 1116 sum_temp 1117 1118 WRITE (unit_nr, '(T2,A)') REPEAT("-", 79) 1119 WRITE (unit_nr, '()') 1120 1121 IF (almo_scf_env%ndomains .LE. 64) THEN 1122 1123 ! print fragment info 1124 WRITE (unit_nr, '(T2,A10,A13,A13,A13,A13,A13)') & 1125 "Fragment", "Basis Set", "Occupied", "Virtual", "Charge", "Deloc Neig" !,"Discarded Virt" 1126 WRITE (unit_nr, '(T2,A)') REPEAT("-", 79) 1127 DO idomain = 1, almo_scf_env%ndomains 1128 1129 SELECT CASE (almo_scf_env%deloc_method) 1130 CASE (almo_deloc_none) 1131 neig_string = "NONE" 1132 CASE (almo_deloc_x, almo_deloc_scf, almo_deloc_x_then_scf) 1133 neig_string = "ALL" 1134 CASE (almo_deloc_xalmo_1diag, almo_deloc_xalmo_x, almo_deloc_xalmo_scf) 1135 WRITE (neig_string, '(I13)') nneighbors(idomain) 1136 CASE DEFAULT 1137 neig_string = "N/A" 1138 END SELECT 1139 1140 WRITE (unit_nr, '(T2,I10,I13,I13,I13,I13,A13)') & 1141 idomain, almo_scf_env%nbasis_of_domain(idomain), & 1142 SUM(almo_scf_env%nocc_of_domain(idomain, :)), & 1143 SUM(almo_scf_env%nvirt_of_domain(idomain, :)), & 1144 !SUM(almo_scf_env%nvirt_disc_of_domain(idomain,:)),& 1145 almo_scf_env%charge_of_domain(idomain), & 1146 ADJUSTR(TRIM(neig_string)) 1147 1148 ENDDO ! cycle over domains 1149 1150 SELECT CASE (almo_scf_env%deloc_method) 1151 CASE (almo_deloc_xalmo_1diag, almo_deloc_xalmo_x, almo_deloc_xalmo_scf) 1152 1153 WRITE (unit_nr, '(T2,A)') REPEAT("-", 79) 1154 1155 ! print fragment neighbors 1156 WRITE (unit_nr, '(T2,A78)') & 1157 "Neighbor lists (including self)" 1158 WRITE (unit_nr, '(T2,A)') REPEAT("-", 79) 1159 DO idomain = 1, almo_scf_env%ndomains 1160 1161 IF (idomain .EQ. 1) THEN 1162 index1_prev = 1 1163 ELSE 1164 index1_prev = almo_scf_env%domain_map(1)%index1(idomain - 1) 1165 ENDIF 1166 1167 WRITE (unit_nr, '(T2,I10,":")') idomain 1168 WRITE (unit_nr, '(T12,11I6)') & 1169 almo_scf_env%domain_map(1)%pairs & 1170 (index1_prev:almo_scf_env%domain_map(1)%index1(idomain) - 1, 1) ! includes self 1171 1172 ENDDO ! cycle over domains 1173 1174 END SELECT 1175 1176 ELSE ! too big to print details for each fragment 1177 1178 WRITE (unit_nr, '(T2,A)') "The system is too big to print details for each fragment." 1179 1180 ENDIF ! how many fragments? 1181 1182 WRITE (unit_nr, '(T2,A)') REPEAT("-", 79) 1183 1184 WRITE (unit_nr, '()') 1185 1186 DEALLOCATE (nneighbors) 1187 1188 ENDIF ! unit_nr > 0 1189 1190 CALL timestop(handle) 1191 1192 END SUBROUTINE almo_scf_print_job_info 1193 1194! ************************************************************************************************** 1195!> \brief Initializes the ALMO SCF copy of the AO overlap matrix 1196!> and all necessary functions (sqrt, inverse...) 1197!> \param matrix_s ... 1198!> \param almo_scf_env ... 1199!> \par History 1200!> 2011.06 created [Rustam Z Khaliullin] 1201!> \author Rustam Z Khaliullin 1202! ************************************************************************************************** 1203 SUBROUTINE almo_scf_init_ao_overlap(matrix_s, almo_scf_env) 1204 TYPE(dbcsr_type), INTENT(IN) :: matrix_s 1205 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 1206 1207 CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_init_ao_overlap', & 1208 routineP = moduleN//':'//routineN 1209 1210 INTEGER :: handle, unit_nr 1211 TYPE(cp_logger_type), POINTER :: logger 1212 1213 CALL timeset(routineN, handle) 1214 1215 ! get a useful output_unit 1216 logger => cp_get_default_logger() 1217 IF (logger%para_env%ionode) THEN 1218 unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.) 1219 ELSE 1220 unit_nr = -1 1221 ENDIF 1222 1223 ! make almo copy of S 1224 ! also copy S to S_blk (i.e. to S with the domain structure imposed) 1225 IF (almo_scf_env%orthogonal_basis) THEN 1226 CALL dbcsr_set(almo_scf_env%matrix_s(1), 0.0_dp) 1227 CALL dbcsr_add_on_diag(almo_scf_env%matrix_s(1), 1.0_dp) 1228 CALL dbcsr_set(almo_scf_env%matrix_s_blk(1), 0.0_dp) 1229 CALL dbcsr_add_on_diag(almo_scf_env%matrix_s_blk(1), 1.0_dp) 1230 ELSE 1231 CALL matrix_qs_to_almo(matrix_s, almo_scf_env%matrix_s(1), & 1232 almo_scf_env%mat_distr_aos, .FALSE.) 1233 CALL dbcsr_copy(almo_scf_env%matrix_s_blk(1), & 1234 almo_scf_env%matrix_s(1), keep_sparsity=.TRUE.) 1235 ENDIF 1236 1237 CALL dbcsr_filter(almo_scf_env%matrix_s(1), almo_scf_env%eps_filter) 1238 CALL dbcsr_filter(almo_scf_env%matrix_s_blk(1), almo_scf_env%eps_filter) 1239 1240 IF (almo_scf_env%almo_update_algorithm .EQ. almo_scf_diag) THEN 1241 CALL matrix_sqrt_Newton_Schulz(almo_scf_env%matrix_s_blk_sqrt(1), & 1242 almo_scf_env%matrix_s_blk_sqrt_inv(1), & 1243 almo_scf_env%matrix_s_blk(1), & 1244 threshold=almo_scf_env%eps_filter, & 1245 order=almo_scf_env%order_lanczos, & 1246 !order=0, & 1247 eps_lanczos=almo_scf_env%eps_lanczos, & 1248 max_iter_lanczos=almo_scf_env%max_iter_lanczos) 1249 ELSE IF (almo_scf_env%almo_update_algorithm .EQ. almo_scf_dm_sign) THEN 1250 CALL invert_Hotelling(almo_scf_env%matrix_s_blk_inv(1), & 1251 almo_scf_env%matrix_s_blk(1), & 1252 threshold=almo_scf_env%eps_filter, & 1253 filter_eps=almo_scf_env%eps_filter) 1254 ENDIF 1255 1256 CALL timestop(handle) 1257 1258 END SUBROUTINE almo_scf_init_ao_overlap 1259 1260! ************************************************************************************************** 1261!> \brief Selects the subroutine for the optimization of block-daigonal ALMOs. 1262!> Keep it short and clean. 1263!> \param qs_env ... 1264!> \param almo_scf_env ... 1265!> \par History 1266!> 2011.11 created [Rustam Z Khaliullin] 1267!> \author Rustam Z Khaliullin 1268! ************************************************************************************************** 1269 SUBROUTINE almo_scf_main(qs_env, almo_scf_env) 1270 TYPE(qs_environment_type), POINTER :: qs_env 1271 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 1272 1273 CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_main', routineP = moduleN//':'//routineN 1274 1275 INTEGER :: handle, ispin, unit_nr 1276 TYPE(cp_logger_type), POINTER :: logger 1277 1278 CALL timeset(routineN, handle) 1279 1280 ! get a useful output_unit 1281 logger => cp_get_default_logger() 1282 IF (logger%para_env%ionode) THEN 1283 unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.) 1284 ELSE 1285 unit_nr = -1 1286 ENDIF 1287 1288 SELECT CASE (almo_scf_env%almo_update_algorithm) 1289 CASE (almo_scf_pcg) 1290 1291 ! ALMO PCG optimizer as a special case of XALMO PCG 1292 CALL almo_scf_xalmo_pcg(qs_env=qs_env, & 1293 almo_scf_env=almo_scf_env, & 1294 optimizer=almo_scf_env%opt_block_diag_pcg, & 1295 quench_t=almo_scf_env%quench_t_blk, & 1296 matrix_t_in=almo_scf_env%matrix_t_blk, & 1297 matrix_t_out=almo_scf_env%matrix_t_blk, & 1298 assume_t0_q0x=.FALSE., & 1299 perturbation_only=.FALSE., & 1300 special_case=xalmo_case_block_diag) 1301 1302 DO ispin = 1, almo_scf_env%nspins 1303 CALL orthogonalize_mos(ket=almo_scf_env%matrix_t_blk(ispin), & 1304 overlap=almo_scf_env%matrix_sigma_blk(ispin), & 1305 metric=almo_scf_env%matrix_s_blk(1), & 1306 retain_locality=.TRUE., & 1307 only_normalize=.FALSE., & 1308 nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), & 1309 eps_filter=almo_scf_env%eps_filter, & 1310 order_lanczos=almo_scf_env%order_lanczos, & 1311 eps_lanczos=almo_scf_env%eps_lanczos, & 1312 max_iter_lanczos=almo_scf_env%max_iter_lanczos) 1313 ENDDO 1314 !CALL almo_scf_t_blk_to_t_blk_orthonormal(almo_scf_env) 1315 1316 CASE (almo_scf_diag) 1317 1318 ! mixing/DIIS optimizer 1319 CALL almo_scf_block_diagonal(qs_env, almo_scf_env, & 1320 almo_scf_env%opt_block_diag_diis) 1321 1322 CASE (almo_scf_skip) 1323 1324 DO ispin = 1, almo_scf_env%nspins 1325 CALL orthogonalize_mos(ket=almo_scf_env%matrix_t_blk(ispin), & 1326 overlap=almo_scf_env%matrix_sigma_blk(ispin), & 1327 metric=almo_scf_env%matrix_s_blk(1), & 1328 retain_locality=.TRUE., & 1329 only_normalize=.FALSE., & 1330 nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), & 1331 eps_filter=almo_scf_env%eps_filter, & 1332 order_lanczos=almo_scf_env%order_lanczos, & 1333 eps_lanczos=almo_scf_env%eps_lanczos, & 1334 max_iter_lanczos=almo_scf_env%max_iter_lanczos) 1335 ENDDO 1336 !CALL almo_scf_t_blk_to_t_blk_orthonormal(almo_scf_env) 1337 1338 END SELECT 1339 1340 ! we might need a copy of the converged KS and sigma_inv 1341 DO ispin = 1, almo_scf_env%nspins 1342 CALL dbcsr_copy(almo_scf_env%matrix_ks_0deloc(ispin), & 1343 almo_scf_env%matrix_ks(ispin)) 1344 CALL dbcsr_copy(almo_scf_env%matrix_sigma_inv_0deloc(ispin), & 1345 almo_scf_env%matrix_sigma_inv(ispin)) 1346 ENDDO 1347 1348 CALL timestop(handle) 1349 1350 END SUBROUTINE almo_scf_main 1351 1352! ************************************************************************************************** 1353!> \brief selects various post scf routines 1354!> \param qs_env ... 1355!> \param almo_scf_env ... 1356!> \par History 1357!> 2011.06 created [Rustam Z Khaliullin] 1358!> \author Rustam Z Khaliullin 1359! ************************************************************************************************** 1360 SUBROUTINE almo_scf_delocalization(qs_env, almo_scf_env) 1361 1362 TYPE(qs_environment_type), POINTER :: qs_env 1363 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 1364 1365 CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_delocalization', & 1366 routineP = moduleN//':'//routineN 1367 1368 INTEGER :: col, handle, hold, iblock_col, & 1369 iblock_row, ispin, mynode, & 1370 nblkcols_tot, nblkrows_tot, row, & 1371 unit_nr 1372 LOGICAL :: almo_experimental, tr 1373 REAL(KIND=dp), DIMENSION(:, :), POINTER :: p_new_block 1374 TYPE(cp_logger_type), POINTER :: logger 1375 TYPE(dbcsr_distribution_type) :: dist 1376 TYPE(dbcsr_type), ALLOCATABLE, DIMENSION(:) :: no_quench 1377 TYPE(optimizer_options_type) :: arbitrary_optimizer 1378 1379 CALL timeset(routineN, handle) 1380 1381 ! get a useful output_unit 1382 logger => cp_get_default_logger() 1383 IF (logger%para_env%ionode) THEN 1384 unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.) 1385 ELSE 1386 unit_nr = -1 1387 ENDIF 1388 1389 ! create a local optimizer that handles XALMO DIIS 1390 ! the options of this optimizer are arbitrary because 1391 ! XALMO DIIS SCF does not converge for yet unknown reasons and 1392 ! currently used in the code to get perturbative estimates only 1393 arbitrary_optimizer%optimizer_type = optimizer_diis 1394 arbitrary_optimizer%max_iter = 3 1395 arbitrary_optimizer%eps_error = 1.0E-6_dp 1396 arbitrary_optimizer%ndiis = 2 1397 1398 SELECT CASE (almo_scf_env%deloc_method) 1399 CASE (almo_deloc_x, almo_deloc_scf, almo_deloc_x_then_scf) 1400 1401 ! RZK-warning hack into the quenched routine: 1402 ! create a quench matrix with all-all-all blocks 1.0 1403 ! it is a waste of memory but since matrices are distributed 1404 ! we can tolerate it for now 1405 ALLOCATE (no_quench(almo_scf_env%nspins)) 1406 CALL dbcsr_create(no_quench(1), & 1407 template=almo_scf_env%matrix_t(1), & 1408 matrix_type=dbcsr_type_no_symmetry) 1409 CALL dbcsr_get_info(no_quench(1), distribution=dist) 1410 CALL dbcsr_distribution_get(dist, mynode=mynode) 1411 CALL dbcsr_work_create(no_quench(1), & 1412 work_mutable=.TRUE.) 1413 nblkrows_tot = dbcsr_nblkrows_total(no_quench(1)) 1414 nblkcols_tot = dbcsr_nblkcols_total(no_quench(1)) 1415 ! RZK-warning: is it a quadratic-scaling routine? 1416 ! As a matter of fact it is! But this block treats 1417 ! fully delocalized MOs. So it is unavoidable. 1418 ! C'est la vie 1419 DO row = 1, nblkrows_tot 1420 DO col = 1, nblkcols_tot 1421 tr = .FALSE. 1422 iblock_row = row 1423 iblock_col = col 1424 CALL dbcsr_get_stored_coordinates(no_quench(1), & 1425 iblock_row, iblock_col, hold) 1426 IF (hold .EQ. mynode) THEN 1427 NULLIFY (p_new_block) 1428 CALL dbcsr_reserve_block2d(no_quench(1), & 1429 iblock_row, iblock_col, p_new_block) 1430 CPASSERT(ASSOCIATED(p_new_block)) 1431 p_new_block(:, :) = 1.0_dp 1432 ENDIF 1433 ENDDO 1434 ENDDO 1435 CALL dbcsr_finalize(no_quench(1)) 1436 IF (almo_scf_env%nspins .GT. 1) THEN 1437 DO ispin = 2, almo_scf_env%nspins 1438 CALL dbcsr_create(no_quench(ispin), & 1439 template=almo_scf_env%matrix_t(1), & 1440 matrix_type=dbcsr_type_no_symmetry) 1441 CALL dbcsr_copy(no_quench(ispin), no_quench(1)) 1442 ENDDO 1443 ENDIF 1444 1445 END SELECT 1446 1447 SELECT CASE (almo_scf_env%deloc_method) 1448 CASE (almo_deloc_none, almo_deloc_scf) 1449 1450 DO ispin = 1, almo_scf_env%nspins 1451 CALL dbcsr_copy(almo_scf_env%matrix_t(ispin), & 1452 almo_scf_env%matrix_t_blk(ispin)) 1453 ENDDO 1454 1455 CASE (almo_deloc_x, almo_deloc_xk, almo_deloc_x_then_scf) 1456 1457 !!!! RZK-warning a whole class of delocalization methods 1458 !!!! are commented out at the moment because some of their 1459 !!!! routines have not been thoroughly tested. 1460 1461 !!!! if we have virtuals pre-optimize and truncate them 1462 !!!IF (almo_scf_env%need_virtuals) THEN 1463 !!! SELECT CASE (almo_scf_env%deloc_truncate_virt) 1464 !!! CASE (virt_full) 1465 !!! ! simply copy virtual orbitals from matrix_v_full_blk to matrix_v_blk 1466 !!! DO ispin=1,almo_scf_env%nspins 1467 !!! CALL dbcsr_copy(almo_scf_env%matrix_v_blk(ispin),& 1468 !!! almo_scf_env%matrix_v_full_blk(ispin)) 1469 !!! ENDDO 1470 !!! CASE (virt_number,virt_occ_size) 1471 !!! CALL split_v_blk(almo_scf_env) 1472 !!! !CALL truncate_subspace_v_blk(qs_env,almo_scf_env) 1473 !!! CASE DEFAULT 1474 !!! CPErrorMessage(cp_failure_level,routineP,"illegal method for virtual space truncation") 1475 !!! CPPrecondition(.FALSE.,cp_failure_level,routineP,failure) 1476 !!! END SELECT 1477 !!!ENDIF 1478 !!!CALL harris_foulkes_correction(qs_env,almo_scf_env) 1479 1480 IF (almo_scf_env%xalmo_update_algorithm .EQ. almo_scf_pcg) THEN 1481 1482 CALL almo_scf_xalmo_pcg(qs_env=qs_env, & 1483 almo_scf_env=almo_scf_env, & 1484 optimizer=almo_scf_env%opt_xalmo_pcg, & 1485 quench_t=no_quench, & 1486 matrix_t_in=almo_scf_env%matrix_t_blk, & 1487 matrix_t_out=almo_scf_env%matrix_t, & 1488 assume_t0_q0x=(almo_scf_env%xalmo_trial_wf .EQ. xalmo_trial_r0_out), & 1489 perturbation_only=.TRUE., & 1490 special_case=xalmo_case_fully_deloc) 1491 1492 ELSE IF (almo_scf_env%xalmo_update_algorithm .EQ. almo_scf_trustr) THEN 1493 1494 CALL almo_scf_xalmo_trustr(qs_env=qs_env, & 1495 almo_scf_env=almo_scf_env, & 1496 optimizer=almo_scf_env%opt_xalmo_trustr, & 1497 quench_t=no_quench, & 1498 matrix_t_in=almo_scf_env%matrix_t_blk, & 1499 matrix_t_out=almo_scf_env%matrix_t, & 1500 perturbation_only=.TRUE., & 1501 special_case=xalmo_case_fully_deloc) 1502 1503 ELSE 1504 1505 CPABORT("Other algorithms do not exist") 1506 1507 ENDIF 1508 1509 CASE (almo_deloc_xalmo_1diag) 1510 1511 IF (almo_scf_env%xalmo_update_algorithm .EQ. almo_scf_diag) THEN 1512 1513 almo_scf_env%perturbative_delocalization = .TRUE. 1514 DO ispin = 1, almo_scf_env%nspins 1515 CALL dbcsr_copy(almo_scf_env%matrix_t(ispin), & 1516 almo_scf_env%matrix_t_blk(ispin)) 1517 ENDDO 1518 CALL almo_scf_xalmo_eigensolver(qs_env, almo_scf_env, & 1519 arbitrary_optimizer) 1520 1521 ELSE 1522 1523 CPABORT("Other algorithms do not exist") 1524 1525 ENDIF 1526 1527 CASE (almo_deloc_xalmo_x) 1528 1529 IF (almo_scf_env%xalmo_update_algorithm .EQ. almo_scf_pcg) THEN 1530 1531 CALL almo_scf_xalmo_pcg(qs_env=qs_env, & 1532 almo_scf_env=almo_scf_env, & 1533 optimizer=almo_scf_env%opt_xalmo_pcg, & 1534 quench_t=almo_scf_env%quench_t, & 1535 matrix_t_in=almo_scf_env%matrix_t_blk, & 1536 matrix_t_out=almo_scf_env%matrix_t, & 1537 assume_t0_q0x=(almo_scf_env%xalmo_trial_wf .EQ. xalmo_trial_r0_out), & 1538 perturbation_only=.TRUE., & 1539 special_case=xalmo_case_normal) 1540 1541 ELSE IF (almo_scf_env%xalmo_update_algorithm .EQ. almo_scf_trustr) THEN 1542 1543 CALL almo_scf_xalmo_trustr(qs_env=qs_env, & 1544 almo_scf_env=almo_scf_env, & 1545 optimizer=almo_scf_env%opt_xalmo_trustr, & 1546 quench_t=almo_scf_env%quench_t, & 1547 matrix_t_in=almo_scf_env%matrix_t_blk, & 1548 matrix_t_out=almo_scf_env%matrix_t, & 1549 perturbation_only=.TRUE., & 1550 special_case=xalmo_case_normal) 1551 1552 ELSE 1553 1554 CPABORT("Other algorithms do not exist") 1555 1556 ENDIF 1557 1558 CASE (almo_deloc_xalmo_scf) 1559 1560 IF (almo_scf_env%xalmo_update_algorithm .EQ. almo_scf_diag) THEN 1561 1562 CPABORT("Should not be here: convergence will fail!") 1563 1564 almo_scf_env%perturbative_delocalization = .FALSE. 1565 DO ispin = 1, almo_scf_env%nspins 1566 CALL dbcsr_copy(almo_scf_env%matrix_t(ispin), & 1567 almo_scf_env%matrix_t_blk(ispin)) 1568 ENDDO 1569 CALL almo_scf_xalmo_eigensolver(qs_env, almo_scf_env, & 1570 arbitrary_optimizer) 1571 1572 ELSE IF (almo_scf_env%xalmo_update_algorithm .EQ. almo_scf_pcg) THEN 1573 1574 CALL almo_scf_xalmo_pcg(qs_env=qs_env, & 1575 almo_scf_env=almo_scf_env, & 1576 optimizer=almo_scf_env%opt_xalmo_pcg, & 1577 quench_t=almo_scf_env%quench_t, & 1578 matrix_t_in=almo_scf_env%matrix_t_blk, & 1579 matrix_t_out=almo_scf_env%matrix_t, & 1580 assume_t0_q0x=(almo_scf_env%xalmo_trial_wf .EQ. xalmo_trial_r0_out), & 1581 perturbation_only=.FALSE., & 1582 special_case=xalmo_case_normal) 1583 1584 ! RZK-warning THIS IS A HACK TO GET ORBITAL ENERGIES 1585 almo_experimental = .FALSE. 1586 IF (almo_experimental) THEN 1587 almo_scf_env%perturbative_delocalization = .TRUE. 1588 !DO ispin=1,almo_scf_env%nspins 1589 ! CALL dbcsr_copy(almo_scf_env%matrix_t(ispin),& 1590 ! almo_scf_env%matrix_t_blk(ispin)) 1591 !ENDDO 1592 CALL almo_scf_xalmo_eigensolver(qs_env, almo_scf_env, & 1593 arbitrary_optimizer) 1594 ENDIF ! experimental 1595 1596 ELSE IF (almo_scf_env%xalmo_update_algorithm .EQ. almo_scf_trustr) THEN 1597 1598 CALL almo_scf_xalmo_trustr(qs_env=qs_env, & 1599 almo_scf_env=almo_scf_env, & 1600 optimizer=almo_scf_env%opt_xalmo_trustr, & 1601 quench_t=almo_scf_env%quench_t, & 1602 matrix_t_in=almo_scf_env%matrix_t_blk, & 1603 matrix_t_out=almo_scf_env%matrix_t, & 1604 perturbation_only=.FALSE., & 1605 special_case=xalmo_case_normal) 1606 1607 ELSE 1608 1609 CPABORT("Other algorithms do not exist") 1610 1611 ENDIF 1612 1613 CASE DEFAULT 1614 1615 CPABORT("Illegal delocalization method") 1616 1617 END SELECT 1618 1619 SELECT CASE (almo_scf_env%deloc_method) 1620 CASE (almo_deloc_scf, almo_deloc_x_then_scf) 1621 1622 IF (almo_scf_env%deloc_truncate_virt .NE. virt_full) THEN 1623 CPABORT("full scf is NYI for truncated virtual space") 1624 ENDIF 1625 1626 IF (almo_scf_env%xalmo_update_algorithm .EQ. almo_scf_pcg) THEN 1627 1628 CALL almo_scf_xalmo_pcg(qs_env=qs_env, & 1629 almo_scf_env=almo_scf_env, & 1630 optimizer=almo_scf_env%opt_xalmo_pcg, & 1631 quench_t=no_quench, & 1632 matrix_t_in=almo_scf_env%matrix_t, & 1633 matrix_t_out=almo_scf_env%matrix_t, & 1634 assume_t0_q0x=.FALSE., & 1635 perturbation_only=.FALSE., & 1636 special_case=xalmo_case_fully_deloc) 1637 1638 ELSE IF (almo_scf_env%xalmo_update_algorithm .EQ. almo_scf_trustr) THEN 1639 1640 CALL almo_scf_xalmo_trustr(qs_env=qs_env, & 1641 almo_scf_env=almo_scf_env, & 1642 optimizer=almo_scf_env%opt_xalmo_trustr, & 1643 quench_t=no_quench, & 1644 matrix_t_in=almo_scf_env%matrix_t, & 1645 matrix_t_out=almo_scf_env%matrix_t, & 1646 perturbation_only=.FALSE., & 1647 special_case=xalmo_case_fully_deloc) 1648 1649 ELSE 1650 1651 CPABORT("Other algorithms do not exist") 1652 1653 ENDIF 1654 1655 END SELECT 1656 1657 ! clean up 1658 SELECT CASE (almo_scf_env%deloc_method) 1659 CASE (almo_deloc_x, almo_deloc_scf, almo_deloc_x_then_scf) 1660 DO ispin = 1, almo_scf_env%nspins 1661 CALL dbcsr_release(no_quench(ispin)) 1662 ENDDO 1663 DEALLOCATE (no_quench) 1664 END SELECT 1665 1666 CALL timestop(handle) 1667 1668 END SUBROUTINE almo_scf_delocalization 1669 1670! ************************************************************************************************** 1671!> \brief orbital localization 1672!> \param qs_env ... 1673!> \param almo_scf_env ... 1674!> \par History 1675!> 2018.09 created [Ziling Luo] 1676!> \author Ziling Luo 1677! ************************************************************************************************** 1678 SUBROUTINE construct_nlmos(qs_env, almo_scf_env) 1679 1680 TYPE(qs_environment_type), POINTER :: qs_env 1681 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 1682 1683 CHARACTER(len=*), PARAMETER :: routineN = 'construct_nlmos', & 1684 routineP = moduleN//':'//routineN 1685 1686 INTEGER :: ispin 1687 1688 IF (almo_scf_env%construct_nlmos) THEN 1689 1690 DO ispin = 1, almo_scf_env%nspins 1691 1692 CALL orthogonalize_mos(ket=almo_scf_env%matrix_t(ispin), & 1693 overlap=almo_scf_env%matrix_sigma(ispin), & 1694 metric=almo_scf_env%matrix_s(1), & 1695 retain_locality=.FALSE., & 1696 only_normalize=.FALSE., & 1697 nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), & 1698 eps_filter=almo_scf_env%eps_filter, & 1699 order_lanczos=almo_scf_env%order_lanczos, & 1700 eps_lanczos=almo_scf_env%eps_lanczos, & 1701 max_iter_lanczos=almo_scf_env%max_iter_lanczos) 1702 ENDDO 1703 1704 CALL construct_nlmos_wrapper(qs_env, almo_scf_env, virtuals=.FALSE.) 1705 1706 IF (almo_scf_env%opt_nlmo_pcg%opt_penalty%virtual_nlmos) THEN 1707 CALL construct_virtuals(almo_scf_env) 1708 CALL construct_nlmos_wrapper(qs_env, almo_scf_env, virtuals=.TRUE.) 1709 ENDIF 1710 1711 IF (almo_scf_env%opt_nlmo_pcg%opt_penalty%compactification_filter_start .GT. 0.0_dp) THEN 1712 CALL nlmo_compactification(qs_env, almo_scf_env, almo_scf_env%matrix_t) 1713 ENDIF 1714 1715 ENDIF 1716 1717 END SUBROUTINE construct_nlmos 1718 1719! ************************************************************************************************** 1720!> \brief Calls NLMO optimization 1721!> \param qs_env ... 1722!> \param almo_scf_env ... 1723!> \param virtuals ... 1724!> \par History 1725!> 2019.10 created [Ziling Luo] 1726!> \author Ziling Luo 1727! ************************************************************************************************** 1728 SUBROUTINE construct_nlmos_wrapper(qs_env, almo_scf_env, virtuals) 1729 1730 TYPE(qs_environment_type), POINTER :: qs_env 1731 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 1732 LOGICAL, INTENT(IN) :: virtuals 1733 1734 CHARACTER(len=*), PARAMETER :: routineN = 'construct_nlmos_wrapper', & 1735 routineP = moduleN//':'//routineN 1736 1737 REAL(KIND=dp) :: det_diff, prev_determinant 1738 1739 almo_scf_env%overlap_determinant = 1.0 1740 ! KEEP: initial_vol_coeff = almo_scf_env%opt_nlmo_pcg%opt_penalty%penalty_strength 1741 almo_scf_env%opt_nlmo_pcg%opt_penalty%penalty_strength = & 1742 -1.0_dp*almo_scf_env%opt_nlmo_pcg%opt_penalty%penalty_strength !NEW1 1743 1744 ! loop over the strength of the orthogonalization penalty 1745 prev_determinant = 10.0_dp 1746 DO WHILE (almo_scf_env%overlap_determinant .GT. almo_scf_env%opt_nlmo_pcg%opt_penalty%final_determinant) 1747 1748 IF (.NOT. virtuals) THEN 1749 CALL almo_scf_construct_nlmos(qs_env=qs_env, & 1750 optimizer=almo_scf_env%opt_nlmo_pcg, & 1751 matrix_s=almo_scf_env%matrix_s(1), & 1752 matrix_mo_in=almo_scf_env%matrix_t, & 1753 matrix_mo_out=almo_scf_env%matrix_t, & 1754 template_matrix_sigma=almo_scf_env%matrix_sigma_inv, & 1755 overlap_determinant=almo_scf_env%overlap_determinant, & 1756 mat_distr_aos=almo_scf_env%mat_distr_aos, & 1757 virtuals=virtuals, & 1758 eps_filter=almo_scf_env%eps_filter) 1759 ELSE 1760 CALL almo_scf_construct_nlmos(qs_env=qs_env, & 1761 optimizer=almo_scf_env%opt_nlmo_pcg, & 1762 matrix_s=almo_scf_env%matrix_s(1), & 1763 matrix_mo_in=almo_scf_env%matrix_v, & 1764 matrix_mo_out=almo_scf_env%matrix_v, & 1765 template_matrix_sigma=almo_scf_env%matrix_sigma_vv, & 1766 overlap_determinant=almo_scf_env%overlap_determinant, & 1767 mat_distr_aos=almo_scf_env%mat_distr_aos, & 1768 virtuals=virtuals, & 1769 eps_filter=almo_scf_env%eps_filter) 1770 1771 ENDIF 1772 1773 det_diff = prev_determinant - almo_scf_env%overlap_determinant 1774 almo_scf_env%opt_nlmo_pcg%opt_penalty%penalty_strength = & 1775 almo_scf_env%opt_nlmo_pcg%opt_penalty%penalty_strength/ & 1776 ABS(almo_scf_env%opt_nlmo_pcg%opt_penalty%penalty_strength_dec_factor) 1777 1778 IF (det_diff < almo_scf_env%opt_nlmo_pcg%opt_penalty%determinant_tolerance) THEN 1779 EXIT 1780 ENDIF 1781 prev_determinant = almo_scf_env%overlap_determinant 1782 1783 ENDDO 1784 1785 END SUBROUTINE construct_nlmos_wrapper 1786 1787! ************************************************************************************************** 1788!> \brief Construct virtual orbitals 1789!> \param almo_scf_env ... 1790!> \par History 1791!> 2019.10 created [Ziling Luo] 1792!> \author Ziling Luo 1793! ************************************************************************************************** 1794 SUBROUTINE construct_virtuals(almo_scf_env) 1795 1796 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 1797 1798 CHARACTER(len=*), PARAMETER :: routineN = 'construct_virtuals', & 1799 routineP = moduleN//':'//routineN 1800 1801 INTEGER :: ispin, n 1802 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: eigenvalues 1803 TYPE(dbcsr_type) :: tempNV1, tempVOcc1, tempVOcc2, tempVV1, & 1804 tempVV2 1805 1806 DO ispin = 1, almo_scf_env%nspins 1807 1808 CALL dbcsr_create(tempNV1, & 1809 template=almo_scf_env%matrix_v(ispin), & 1810 matrix_type=dbcsr_type_no_symmetry) 1811 CALL dbcsr_create(tempVOcc1, & 1812 template=almo_scf_env%matrix_vo(ispin), & 1813 matrix_type=dbcsr_type_no_symmetry) 1814 CALL dbcsr_create(tempVOcc2, & 1815 template=almo_scf_env%matrix_vo(ispin), & 1816 matrix_type=dbcsr_type_no_symmetry) 1817 CALL dbcsr_create(tempVV1, & 1818 template=almo_scf_env%matrix_sigma_vv(ispin), & 1819 matrix_type=dbcsr_type_no_symmetry) 1820 CALL dbcsr_create(tempVV2, & 1821 template=almo_scf_env%matrix_sigma_vv(ispin), & 1822 matrix_type=dbcsr_type_no_symmetry) 1823 1824 ! Generate random virtual matrix 1825 CALL dbcsr_init_random(almo_scf_env%matrix_v(ispin), & 1826 keep_sparsity=.FALSE.) 1827 1828 ! Project the orbital subspace out 1829 CALL dbcsr_multiply("N", "N", 1.0_dp, & 1830 almo_scf_env%matrix_s(1), & 1831 almo_scf_env%matrix_v(ispin), & 1832 0.0_dp, tempNV1, & 1833 filter_eps=almo_scf_env%eps_filter) 1834 1835 CALL dbcsr_multiply("T", "N", 1.0_dp, & 1836 tempNV1, & 1837 almo_scf_env%matrix_t(ispin), & 1838 0.0_dp, tempVOcc1, & 1839 filter_eps=almo_scf_env%eps_filter) 1840 1841 CALL dbcsr_multiply("N", "N", 1.0_dp, & 1842 tempVOcc1, & 1843 almo_scf_env%matrix_sigma_inv(ispin), & 1844 0.0_dp, tempVOcc2, & 1845 filter_eps=almo_scf_env%eps_filter) 1846 1847 CALL dbcsr_multiply("N", "T", 1.0_dp, & 1848 almo_scf_env%matrix_t(ispin), & 1849 tempVOcc2, & 1850 0.0_dp, tempNV1, & 1851 filter_eps=almo_scf_env%eps_filter) 1852 1853 CALL dbcsr_add(almo_scf_env%matrix_v(ispin), tempNV1, 1.0_dp, -1.0_dp) 1854 1855 ! compute VxV overlap 1856 CALL dbcsr_multiply("N", "N", 1.0_dp, & 1857 almo_scf_env%matrix_s(1), & 1858 almo_scf_env%matrix_v(ispin), & 1859 0.0_dp, tempNV1, & 1860 filter_eps=almo_scf_env%eps_filter) 1861 1862 CALL dbcsr_multiply("T", "N", 1.0_dp, & 1863 almo_scf_env%matrix_v(ispin), & 1864 tempNV1, & 1865 0.0_dp, tempVV1, & 1866 filter_eps=almo_scf_env%eps_filter) 1867 1868 CALL orthogonalize_mos(ket=almo_scf_env%matrix_v(ispin), & 1869 overlap=tempVV1, & 1870 metric=almo_scf_env%matrix_s(1), & 1871 retain_locality=.FALSE., & 1872 only_normalize=.FALSE., & 1873 nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), & 1874 eps_filter=almo_scf_env%eps_filter, & 1875 order_lanczos=almo_scf_env%order_lanczos, & 1876 eps_lanczos=almo_scf_env%eps_lanczos, & 1877 max_iter_lanczos=almo_scf_env%max_iter_lanczos) 1878 1879 ! compute VxV block of the KS matrix 1880 CALL dbcsr_multiply("N", "N", 1.0_dp, & 1881 almo_scf_env%matrix_ks(ispin), & 1882 almo_scf_env%matrix_v(ispin), & 1883 0.0_dp, tempNV1, & 1884 filter_eps=almo_scf_env%eps_filter) 1885 1886 CALL dbcsr_multiply("T", "N", 1.0_dp, & 1887 almo_scf_env%matrix_v(ispin), & 1888 tempNV1, & 1889 0.0_dp, tempVV1, & 1890 filter_eps=almo_scf_env%eps_filter) 1891 1892 CALL dbcsr_get_info(tempVV1, nfullrows_total=n) 1893 ALLOCATE (eigenvalues(n)) 1894 CALL cp_dbcsr_syevd(tempVV1, tempVV2, & 1895 eigenvalues, & 1896 para_env=almo_scf_env%para_env, & 1897 blacs_env=almo_scf_env%blacs_env) 1898 DEALLOCATE (eigenvalues) 1899 1900 CALL dbcsr_multiply("N", "N", 1.0_dp, & 1901 almo_scf_env%matrix_v(ispin), & 1902 tempVV2, & 1903 0.0_dp, tempNV1, & 1904 filter_eps=almo_scf_env%eps_filter) 1905 1906 CALL dbcsr_copy(almo_scf_env%matrix_v(ispin), tempNV1) 1907 1908 CALL dbcsr_release(tempNV1) 1909 CALL dbcsr_release(tempVOcc1) 1910 CALL dbcsr_release(tempVOcc2) 1911 CALL dbcsr_release(tempVV1) 1912 CALL dbcsr_release(tempVV2) 1913 1914 ENDDO 1915 1916 END SUBROUTINE construct_virtuals 1917 1918! ************************************************************************************************** 1919!> \brief Compactify (set small blocks to zero) orbitals 1920!> \param qs_env ... 1921!> \param almo_scf_env ... 1922!> \param matrix ... 1923!> \par History 1924!> 2019.10 created [Ziling Luo] 1925!> \author Ziling Luo 1926! ************************************************************************************************** 1927 SUBROUTINE nlmo_compactification(qs_env, almo_scf_env, matrix) 1928 1929 TYPE(qs_environment_type), POINTER :: qs_env 1930 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 1931 TYPE(dbcsr_type), ALLOCATABLE, DIMENSION(:), & 1932 INTENT(IN) :: matrix 1933 1934 CHARACTER(len=*), PARAMETER :: routineN = 'nlmo_compactification', & 1935 routineP = moduleN//':'//routineN 1936 1937 INTEGER :: iblock_col, iblock_col_size, iblock_row, & 1938 iblock_row_size, icol, irow, ispin, & 1939 Ncols, Nrows, nspins, para_group, & 1940 unit_nr 1941 LOGICAL :: element_by_element 1942 REAL(KIND=dp) :: energy, eps_local, eps_start, & 1943 max_element, spin_factor 1944 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: occ, retained 1945 REAL(kind=dp), DIMENSION(:, :), POINTER :: data_p 1946 TYPE(cp_logger_type), POINTER :: logger 1947 TYPE(dbcsr_iterator_type) :: iter 1948 TYPE(dbcsr_type), ALLOCATABLE, DIMENSION(:) :: matrix_p_tmp, matrix_t_tmp 1949 1950 ! define the output_unit 1951 logger => cp_get_default_logger() 1952 IF (logger%para_env%ionode) THEN 1953 unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.) 1954 ELSE 1955 unit_nr = -1 1956 ENDIF 1957 1958 nspins = SIZE(matrix) 1959 element_by_element = .FALSE. 1960 1961 IF (nspins .EQ. 1) THEN 1962 spin_factor = 2.0_dp 1963 ELSE 1964 spin_factor = 1.0_dp 1965 ENDIF 1966 1967 ALLOCATE (matrix_t_tmp(nspins)) 1968 ALLOCATE (matrix_p_tmp(nspins)) 1969 ALLOCATE (retained(nspins)) 1970 ALLOCATE (occ(2)) 1971 1972 DO ispin = 1, nspins 1973 1974 ! init temporary storage 1975 CALL dbcsr_create(matrix_t_tmp(ispin), & 1976 template=matrix(ispin), & 1977 matrix_type=dbcsr_type_no_symmetry) 1978 CALL dbcsr_copy(matrix_t_tmp(ispin), matrix(ispin)) 1979 1980 CALL dbcsr_create(matrix_p_tmp(ispin), & 1981 template=almo_scf_env%matrix_p(ispin), & 1982 matrix_type=dbcsr_type_no_symmetry) 1983 CALL dbcsr_copy(matrix_p_tmp(ispin), almo_scf_env%matrix_p(ispin)) 1984 1985 ENDDO 1986 1987 IF (unit_nr > 0) THEN 1988 WRITE (unit_nr, *) 1989 WRITE (unit_nr, '(T2,A)') & 1990 "Energy dependence on the (block-by-block) filtering of the NLMO coefficients" 1991 IF (unit_nr > 0) WRITE (unit_nr, '(T2,A13,A20,A20,A25)') & 1992 "EPS filter", "Occupation Alpha", "Occupation Beta", "Energy" 1993 ENDIF 1994 1995 eps_start = almo_scf_env%opt_nlmo_pcg%opt_penalty%compactification_filter_start 1996 eps_local = MAX(eps_start, 10E-14_dp) 1997 1998 DO 1999 2000 IF (eps_local > 0.11_dp) EXIT 2001 2002 DO ispin = 1, nspins 2003 2004 retained(ispin) = 0 2005 CALL dbcsr_work_create(matrix_t_tmp(ispin), work_mutable=.TRUE.) 2006 CALL dbcsr_iterator_start(iter, matrix_t_tmp(ispin)) 2007 DO WHILE (dbcsr_iterator_blocks_left(iter)) 2008 CALL dbcsr_iterator_next_block(iter, iblock_row, iblock_col, data_p, & 2009 row_size=iblock_row_size, col_size=iblock_col_size) 2010 DO icol = 1, iblock_col_size 2011 2012 IF (element_by_element) THEN 2013 2014 DO irow = 1, iblock_row_size 2015 IF (ABS(data_p(irow, icol)) .LT. eps_local) THEN 2016 data_p(irow, icol) = 0.0_dp 2017 ELSE 2018 retained(ispin) = retained(ispin) + 1 2019 ENDIF 2020 ENDDO 2021 2022 ELSE ! rows are blocked 2023 2024 max_element = 0.0_dp 2025 DO irow = 1, iblock_row_size 2026 IF (ABS(data_p(irow, icol)) .GT. max_element) THEN 2027 max_element = ABS(data_p(irow, icol)) 2028 ENDIF 2029 ENDDO 2030 IF (max_element .LT. eps_local) THEN 2031 DO irow = 1, iblock_row_size 2032 data_p(irow, icol) = 0.0_dp 2033 ENDDO 2034 ELSE 2035 retained(ispin) = retained(ispin) + iblock_row_size 2036 ENDIF 2037 2038 ENDIF ! block rows? 2039 ENDDO ! icol 2040 2041 ENDDO ! iterator 2042 CALL dbcsr_iterator_stop(iter) 2043 CALL dbcsr_finalize(matrix_t_tmp(ispin)) 2044 CALL dbcsr_filter(matrix_t_tmp(ispin), eps_local) 2045 2046 CALL dbcsr_get_info(matrix_t_tmp(ispin), group=para_group, & 2047 nfullrows_total=Nrows, & 2048 nfullcols_total=Ncols) 2049 CALL mp_sum(retained(ispin), para_group) 2050 2051 !devide by the total no. elements 2052 occ(ispin) = retained(ispin)/Nrows/Ncols 2053 2054 ! compute the global projectors (for the density matrix) 2055 CALL almo_scf_t_to_proj( & 2056 t=matrix_t_tmp(ispin), & 2057 p=matrix_p_tmp(ispin), & 2058 eps_filter=almo_scf_env%eps_filter, & 2059 orthog_orbs=.FALSE., & 2060 nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), & 2061 s=almo_scf_env%matrix_s(1), & 2062 sigma=almo_scf_env%matrix_sigma(ispin), & 2063 sigma_inv=almo_scf_env%matrix_sigma_inv(ispin), & 2064 use_guess=.FALSE., & 2065 algorithm=almo_scf_env%sigma_inv_algorithm, & 2066 inv_eps_factor=almo_scf_env%matrix_iter_eps_error_factor, & 2067 inverse_accelerator=almo_scf_env%order_lanczos, & 2068 eps_lanczos=almo_scf_env%eps_lanczos, & 2069 max_iter_lanczos=almo_scf_env%max_iter_lanczos, & 2070 para_env=almo_scf_env%para_env, & 2071 blacs_env=almo_scf_env%blacs_env) 2072 2073 ! compute dm from the projector(s) 2074 CALL dbcsr_scale(matrix_p_tmp(ispin), spin_factor) 2075 2076 ENDDO 2077 2078 ! the KS matrix is updated outside the spin loop 2079 CALL almo_dm_to_almo_ks(qs_env, & 2080 matrix_p_tmp, & 2081 almo_scf_env%matrix_ks, & 2082 energy, & 2083 almo_scf_env%eps_filter, & 2084 almo_scf_env%mat_distr_aos) 2085 2086 IF (nspins .LT. 2) occ(2) = occ(1) 2087 IF (unit_nr > 0) WRITE (unit_nr, '(T2,E13.3,F20.10,F20.10,F25.15)') & 2088 eps_local, occ(1), occ(2), energy 2089 2090 eps_local = 2.0_dp*eps_local 2091 2092 ENDDO 2093 2094 DO ispin = 1, nspins 2095 2096 CALL dbcsr_release(matrix_t_tmp(ispin)) 2097 CALL dbcsr_release(matrix_p_tmp(ispin)) 2098 2099 ENDDO 2100 2101 DEALLOCATE (matrix_t_tmp) 2102 DEALLOCATE (matrix_p_tmp) 2103 DEALLOCATE (occ) 2104 DEALLOCATE (retained) 2105 2106 END SUBROUTINE nlmo_compactification 2107 2108! ***************************************************************************** 2109!> \brief after SCF we have the final density and KS matrices compute various 2110!> post-scf quantities 2111!> \param qs_env ... 2112!> \param almo_scf_env ... 2113!> \par History 2114!> 2015.03 created [Rustam Z. Khaliullin] 2115!> \author Rustam Z. Khaliullin 2116! ************************************************************************************************** 2117 SUBROUTINE almo_scf_post(qs_env, almo_scf_env) 2118 TYPE(qs_environment_type), POINTER :: qs_env 2119 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 2120 2121 CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_post', routineP = moduleN//':'//routineN 2122 2123 INTEGER :: handle, ispin 2124 TYPE(cp_fm_type), POINTER :: mo_coeff 2125 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_w 2126 TYPE(dbcsr_type), ALLOCATABLE, DIMENSION(:) :: matrix_t_processed 2127 TYPE(mo_set_p_type), DIMENSION(:), POINTER :: mos 2128 TYPE(qs_scf_env_type), POINTER :: scf_env 2129 2130 CALL timeset(routineN, handle) 2131 2132 ! store matrices to speed up the next scf run 2133 CALL almo_scf_store_extrapolation_data(almo_scf_env) 2134 2135 ! orthogonalize orbitals before returning them to QS 2136 ALLOCATE (matrix_t_processed(almo_scf_env%nspins)) 2137 !ALLOCATE (matrix_v_processed(almo_scf_env%nspins)) 2138 2139 DO ispin = 1, almo_scf_env%nspins 2140 2141 CALL dbcsr_create(matrix_t_processed(ispin), & 2142 template=almo_scf_env%matrix_t(ispin), & 2143 matrix_type=dbcsr_type_no_symmetry) 2144 2145 CALL dbcsr_copy(matrix_t_processed(ispin), & 2146 almo_scf_env%matrix_t(ispin)) 2147 2148 !CALL dbcsr_create(matrix_v_processed(ispin), & 2149 ! template=almo_scf_env%matrix_v(ispin), & 2150 ! matrix_type=dbcsr_type_no_symmetry) 2151 2152 !CALL dbcsr_copy(matrix_v_processed(ispin), & 2153 ! almo_scf_env%matrix_v(ispin)) 2154 2155 IF (almo_scf_env%return_orthogonalized_mos) THEN 2156 2157 CALL orthogonalize_mos(ket=matrix_t_processed(ispin), & 2158 overlap=almo_scf_env%matrix_sigma(ispin), & 2159 metric=almo_scf_env%matrix_s(1), & 2160 retain_locality=.FALSE., & 2161 only_normalize=.FALSE., & 2162 nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), & 2163 eps_filter=almo_scf_env%eps_filter, & 2164 order_lanczos=almo_scf_env%order_lanczos, & 2165 eps_lanczos=almo_scf_env%eps_lanczos, & 2166 max_iter_lanczos=almo_scf_env%max_iter_lanczos, & 2167 smear=almo_scf_env%smear) 2168 ENDIF 2169 2170 ENDDO 2171 2172 !! RS-WARNING: If smearing ALMO is requested, rescaled fully-occupied orbitals are returned to QS 2173 !! RS-WARNING: Beware that QS will not be informed about electronic entropy. 2174 !! If you want a quick and dirty transfer to QS energy, uncomment the following hack: 2175 !! IF (almo_scf_env%smear) THEN 2176 !! qs_env%energy%kTS = 0.0_dp 2177 !! DO ispin = 1, almo_scf_env%nspins 2178 !! qs_env%energy%kTS = qs_env%energy%kTS + almo_scf_env%kTS(ispin) 2179 !! END DO 2180 !! END IF 2181 2182 ! return orbitals to QS 2183 NULLIFY (mos, mo_coeff, scf_env) 2184 2185 CALL get_qs_env(qs_env, mos=mos, scf_env=scf_env) 2186 2187 DO ispin = 1, almo_scf_env%nspins 2188 ! Currently only fm version of mo_set is usable. 2189 ! First transform the matrix_t to fm version 2190 CALL get_mo_set(mos(ispin)%mo_set, mo_coeff=mo_coeff) 2191 CALL copy_dbcsr_to_fm(matrix_t_processed(ispin), mo_coeff) 2192 CALL dbcsr_release(matrix_t_processed(ispin)) 2193 ENDDO 2194 DO ispin = 1, almo_scf_env%nspins 2195 CALL dbcsr_release(matrix_t_processed(ispin)) 2196 ENDDO 2197 DEALLOCATE (matrix_t_processed) 2198 2199 ! calculate post scf properties 2200 ! CALL almo_post_scf_compute_properties(qs_env, almo_scf_env) 2201 CALL almo_post_scf_compute_properties(qs_env) 2202 2203 ! In normal scf procedure this is done by scf, 2204 ! but it's not called by ALMO 2205 ! best to release it somewhere else 2206 IF (ASSOCIATED(scf_env)) THEN 2207 CALL scf_env_release(scf_env) 2208 ENDIF 2209 2210 ! compute the W matrix if associated 2211 IF (almo_scf_env%calc_forces) THEN 2212 CALL get_qs_env(qs_env, matrix_w=matrix_w) 2213 IF (ASSOCIATED(matrix_w)) THEN 2214 CALL calculate_w_matrix_almo(matrix_w, almo_scf_env) 2215 ELSE 2216 CPABORT("Matrix W is needed but not associated") 2217 ENDIF 2218 ENDIF 2219 2220 CALL timestop(handle) 2221 2222 END SUBROUTINE almo_scf_post 2223 2224! ************************************************************************************************** 2225!> \brief create various matrices 2226!> \param almo_scf_env ... 2227!> \param matrix_s0 ... 2228!> \par History 2229!> 2011.07 created [Rustam Z Khaliullin] 2230!> \author Rustam Z Khaliullin 2231! ************************************************************************************************** 2232 SUBROUTINE almo_scf_env_create_matrices(almo_scf_env, matrix_s0) 2233 2234 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 2235 TYPE(dbcsr_type), INTENT(IN) :: matrix_s0 2236 2237 CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_env_create_matrices', & 2238 routineP = moduleN//':'//routineN 2239 2240 INTEGER :: handle, ispin, nspins 2241 2242 CALL timeset(routineN, handle) 2243 2244 nspins = almo_scf_env%nspins 2245 2246 ! AO overlap matrix and its various functions 2247 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_s(1), & 2248 matrix_qs=matrix_s0, & 2249 almo_scf_env=almo_scf_env, & 2250 name_new="S", & 2251 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_aobasis/), & 2252 symmetry_new=dbcsr_type_symmetric, & 2253 spin_key=0, & 2254 init_domains=.FALSE.) 2255 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_s_blk(1), & 2256 matrix_qs=matrix_s0, & 2257 almo_scf_env=almo_scf_env, & 2258 name_new="S_BLK", & 2259 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_aobasis/), & 2260 symmetry_new=dbcsr_type_symmetric, & 2261 spin_key=0, & 2262 init_domains=.TRUE.) 2263 IF (almo_scf_env%almo_update_algorithm .EQ. almo_scf_diag) THEN 2264 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_s_blk_sqrt_inv(1), & 2265 matrix_qs=matrix_s0, & 2266 almo_scf_env=almo_scf_env, & 2267 name_new="S_BLK_SQRT_INV", & 2268 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_aobasis/), & 2269 symmetry_new=dbcsr_type_symmetric, & 2270 spin_key=0, & 2271 init_domains=.TRUE.) 2272 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_s_blk_sqrt(1), & 2273 matrix_qs=matrix_s0, & 2274 almo_scf_env=almo_scf_env, & 2275 name_new="S_BLK_SQRT", & 2276 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_aobasis/), & 2277 symmetry_new=dbcsr_type_symmetric, & 2278 spin_key=0, & 2279 init_domains=.TRUE.) 2280 ELSE IF (almo_scf_env%almo_update_algorithm .EQ. almo_scf_dm_sign) THEN 2281 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_s_blk_inv(1), & 2282 matrix_qs=matrix_s0, & 2283 almo_scf_env=almo_scf_env, & 2284 name_new="S_BLK_INV", & 2285 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_aobasis/), & 2286 symmetry_new=dbcsr_type_symmetric, & 2287 spin_key=0, & 2288 init_domains=.TRUE.) 2289 ENDIF 2290 2291 ! MO coeff matrices and their derivatives 2292 ALLOCATE (almo_scf_env%matrix_t_blk(nspins)) 2293 ALLOCATE (almo_scf_env%quench_t_blk(nspins)) 2294 ALLOCATE (almo_scf_env%matrix_err_blk(nspins)) 2295 ALLOCATE (almo_scf_env%matrix_err_xx(nspins)) 2296 ALLOCATE (almo_scf_env%matrix_sigma(nspins)) 2297 ALLOCATE (almo_scf_env%matrix_sigma_inv(nspins)) 2298 ALLOCATE (almo_scf_env%matrix_sigma_sqrt(nspins)) 2299 ALLOCATE (almo_scf_env%matrix_sigma_sqrt_inv(nspins)) 2300 ALLOCATE (almo_scf_env%matrix_sigma_blk(nspins)) 2301 ALLOCATE (almo_scf_env%matrix_sigma_inv_0deloc(nspins)) 2302 ALLOCATE (almo_scf_env%matrix_t(nspins)) 2303 ALLOCATE (almo_scf_env%matrix_t_tr(nspins)) 2304 DO ispin = 1, nspins 2305 ! create the blocked quencher 2306 CALL matrix_almo_create(matrix_new=almo_scf_env%quench_t_blk(ispin), & 2307 matrix_qs=matrix_s0, & 2308 almo_scf_env=almo_scf_env, & 2309 name_new="Q_BLK", & 2310 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_occ/), & 2311 symmetry_new=dbcsr_type_no_symmetry, & 2312 spin_key=ispin, & 2313 init_domains=.TRUE.) 2314 ! create ALMO coefficient matrix 2315 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_t_blk(ispin), & 2316 matrix_qs=matrix_s0, & 2317 almo_scf_env=almo_scf_env, & 2318 name_new="T_BLK", & 2319 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_occ/), & 2320 symmetry_new=dbcsr_type_no_symmetry, & 2321 spin_key=ispin, & 2322 init_domains=.TRUE.) 2323 ! create the error matrix 2324 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_err_blk(ispin), & 2325 matrix_qs=matrix_s0, & 2326 almo_scf_env=almo_scf_env, & 2327 name_new="ERR_BLK", & 2328 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_aobasis/), & 2329 symmetry_new=dbcsr_type_no_symmetry, & 2330 spin_key=ispin, & 2331 init_domains=.TRUE.) 2332 ! create the error matrix for the quenched ALMOs 2333 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_err_xx(ispin), & 2334 matrix_qs=matrix_s0, & 2335 almo_scf_env=almo_scf_env, & 2336 name_new="ERR_XX", & 2337 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_occ/), & 2338 symmetry_new=dbcsr_type_no_symmetry, & 2339 spin_key=ispin, & 2340 init_domains=.FALSE.) 2341 ! create a matrix with dimensions of a transposed mo coefficient matrix 2342 ! it might be necessary to perform the correction step using cayley 2343 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_t_tr(ispin), & 2344 matrix_qs=matrix_s0, & 2345 almo_scf_env=almo_scf_env, & 2346 name_new="T_TR", & 2347 size_keys=(/almo_mat_dim_occ, almo_mat_dim_aobasis/), & 2348 symmetry_new=dbcsr_type_no_symmetry, & 2349 spin_key=ispin, & 2350 init_domains=.FALSE.) 2351 ! create mo overlap matrix 2352 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_sigma(ispin), & 2353 matrix_qs=matrix_s0, & 2354 almo_scf_env=almo_scf_env, & 2355 name_new="SIG", & 2356 size_keys=(/almo_mat_dim_occ, almo_mat_dim_occ/), & 2357 symmetry_new=dbcsr_type_symmetric, & 2358 spin_key=ispin, & 2359 init_domains=.FALSE.) 2360 ! create blocked mo overlap matrix 2361 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_sigma_blk(ispin), & 2362 matrix_qs=matrix_s0, & 2363 almo_scf_env=almo_scf_env, & 2364 name_new="SIG_BLK", & 2365 size_keys=(/almo_mat_dim_occ, almo_mat_dim_occ/), & 2366 symmetry_new=dbcsr_type_symmetric, & 2367 spin_key=ispin, & 2368 init_domains=.TRUE.) 2369 ! create blocked inverse mo overlap matrix 2370 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_sigma_inv_0deloc(ispin), & 2371 matrix_qs=matrix_s0, & 2372 almo_scf_env=almo_scf_env, & 2373 name_new="SIGINV_BLK", & 2374 size_keys=(/almo_mat_dim_occ, almo_mat_dim_occ/), & 2375 symmetry_new=dbcsr_type_symmetric, & 2376 spin_key=ispin, & 2377 init_domains=.TRUE.) 2378 ! create inverse mo overlap matrix 2379 CALL matrix_almo_create( & 2380 matrix_new=almo_scf_env%matrix_sigma_inv(ispin), & 2381 matrix_qs=matrix_s0, & 2382 almo_scf_env=almo_scf_env, & 2383 name_new="SIGINV", & 2384 size_keys=(/almo_mat_dim_occ, almo_mat_dim_occ/), & 2385 symmetry_new=dbcsr_type_symmetric, & 2386 spin_key=ispin, & 2387 init_domains=.FALSE.) 2388 ! create various templates that will be necessary later 2389 CALL matrix_almo_create( & 2390 matrix_new=almo_scf_env%matrix_t(ispin), & 2391 matrix_qs=matrix_s0, & 2392 almo_scf_env=almo_scf_env, & 2393 name_new="T", & 2394 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_occ/), & 2395 symmetry_new=dbcsr_type_no_symmetry, & 2396 spin_key=ispin, & 2397 init_domains=.FALSE.) 2398 CALL dbcsr_create(almo_scf_env%matrix_sigma_sqrt(ispin), & 2399 template=almo_scf_env%matrix_sigma(ispin), & 2400 matrix_type=dbcsr_type_no_symmetry) 2401 CALL dbcsr_create(almo_scf_env%matrix_sigma_sqrt_inv(ispin), & 2402 template=almo_scf_env%matrix_sigma(ispin), & 2403 matrix_type=dbcsr_type_no_symmetry) 2404 ENDDO 2405 2406 ! create virtual orbitals if necessary 2407 IF (almo_scf_env%need_virtuals) THEN 2408 ALLOCATE (almo_scf_env%matrix_v_blk(nspins)) 2409 ALLOCATE (almo_scf_env%matrix_v_full_blk(nspins)) 2410 ALLOCATE (almo_scf_env%matrix_v(nspins)) 2411 ALLOCATE (almo_scf_env%matrix_vo(nspins)) 2412 ALLOCATE (almo_scf_env%matrix_x(nspins)) 2413 ALLOCATE (almo_scf_env%matrix_ov(nspins)) 2414 ALLOCATE (almo_scf_env%matrix_ov_full(nspins)) 2415 ALLOCATE (almo_scf_env%matrix_sigma_vv(nspins)) 2416 ALLOCATE (almo_scf_env%matrix_sigma_vv_blk(nspins)) 2417 ALLOCATE (almo_scf_env%matrix_sigma_vv_sqrt(nspins)) 2418 ALLOCATE (almo_scf_env%matrix_sigma_vv_sqrt_inv(nspins)) 2419 ALLOCATE (almo_scf_env%matrix_vv_full_blk(nspins)) 2420 2421 IF (almo_scf_env%deloc_truncate_virt .NE. virt_full) THEN 2422 ALLOCATE (almo_scf_env%matrix_k_blk(nspins)) 2423 ALLOCATE (almo_scf_env%matrix_k_blk_ones(nspins)) 2424 ALLOCATE (almo_scf_env%matrix_k_tr(nspins)) 2425 ALLOCATE (almo_scf_env%matrix_v_disc(nspins)) 2426 ALLOCATE (almo_scf_env%matrix_v_disc_blk(nspins)) 2427 ALLOCATE (almo_scf_env%matrix_ov_disc(nspins)) 2428 ALLOCATE (almo_scf_env%matrix_vv_disc_blk(nspins)) 2429 ALLOCATE (almo_scf_env%matrix_vv_disc(nspins)) 2430 ALLOCATE (almo_scf_env%opt_k_t_dd(nspins)) 2431 ALLOCATE (almo_scf_env%opt_k_t_rr(nspins)) 2432 ALLOCATE (almo_scf_env%opt_k_denom(nspins)) 2433 ENDIF 2434 2435 DO ispin = 1, nspins 2436 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_v_full_blk(ispin), & 2437 matrix_qs=matrix_s0, & 2438 almo_scf_env=almo_scf_env, & 2439 name_new="V_FULL_BLK", & 2440 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_virt_full/), & 2441 symmetry_new=dbcsr_type_no_symmetry, & 2442 spin_key=ispin, & 2443 init_domains=.FALSE.) 2444 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_v_blk(ispin), & 2445 matrix_qs=matrix_s0, & 2446 almo_scf_env=almo_scf_env, & 2447 name_new="V_BLK", & 2448 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_virt/), & 2449 symmetry_new=dbcsr_type_no_symmetry, & 2450 spin_key=ispin, & 2451 init_domains=.FALSE.) 2452 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_v(ispin), & 2453 matrix_qs=matrix_s0, & 2454 almo_scf_env=almo_scf_env, & 2455 name_new="V", & 2456 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_virt/), & 2457 symmetry_new=dbcsr_type_no_symmetry, & 2458 spin_key=ispin, & 2459 init_domains=.FALSE.) 2460 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_ov_full(ispin), & 2461 matrix_qs=matrix_s0, & 2462 almo_scf_env=almo_scf_env, & 2463 name_new="OV_FULL", & 2464 size_keys=(/almo_mat_dim_occ, almo_mat_dim_virt_full/), & 2465 symmetry_new=dbcsr_type_no_symmetry, & 2466 spin_key=ispin, & 2467 init_domains=.FALSE.) 2468 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_ov(ispin), & 2469 matrix_qs=matrix_s0, & 2470 almo_scf_env=almo_scf_env, & 2471 name_new="OV", & 2472 size_keys=(/almo_mat_dim_occ, almo_mat_dim_virt/), & 2473 symmetry_new=dbcsr_type_no_symmetry, & 2474 spin_key=ispin, & 2475 init_domains=.FALSE.) 2476 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_vo(ispin), & 2477 matrix_qs=matrix_s0, & 2478 almo_scf_env=almo_scf_env, & 2479 name_new="VO", & 2480 size_keys=(/almo_mat_dim_virt, almo_mat_dim_occ/), & 2481 symmetry_new=dbcsr_type_no_symmetry, & 2482 spin_key=ispin, & 2483 init_domains=.FALSE.) 2484 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_x(ispin), & 2485 matrix_qs=matrix_s0, & 2486 almo_scf_env=almo_scf_env, & 2487 name_new="VO", & 2488 size_keys=(/almo_mat_dim_virt, almo_mat_dim_occ/), & 2489 symmetry_new=dbcsr_type_no_symmetry, & 2490 spin_key=ispin, & 2491 init_domains=.FALSE.) 2492 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_sigma_vv(ispin), & 2493 matrix_qs=matrix_s0, & 2494 almo_scf_env=almo_scf_env, & 2495 name_new="SIG_VV", & 2496 size_keys=(/almo_mat_dim_virt, almo_mat_dim_virt/), & 2497 symmetry_new=dbcsr_type_symmetric, & 2498 spin_key=ispin, & 2499 init_domains=.FALSE.) 2500 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_vv_full_blk(ispin), & 2501 matrix_qs=matrix_s0, & 2502 almo_scf_env=almo_scf_env, & 2503 name_new="VV_FULL_BLK", & 2504 size_keys=(/almo_mat_dim_virt_full, almo_mat_dim_virt_full/), & 2505 symmetry_new=dbcsr_type_no_symmetry, & 2506 spin_key=ispin, & 2507 init_domains=.TRUE.) 2508 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_sigma_vv_blk(ispin), & 2509 matrix_qs=matrix_s0, & 2510 almo_scf_env=almo_scf_env, & 2511 name_new="SIG_VV_BLK", & 2512 size_keys=(/almo_mat_dim_virt, almo_mat_dim_virt/), & 2513 symmetry_new=dbcsr_type_symmetric, & 2514 spin_key=ispin, & 2515 init_domains=.TRUE.) 2516 CALL dbcsr_create(almo_scf_env%matrix_sigma_vv_sqrt(ispin), & 2517 template=almo_scf_env%matrix_sigma_vv(ispin), & 2518 matrix_type=dbcsr_type_no_symmetry) 2519 CALL dbcsr_create(almo_scf_env%matrix_sigma_vv_sqrt_inv(ispin), & 2520 template=almo_scf_env%matrix_sigma_vv(ispin), & 2521 matrix_type=dbcsr_type_no_symmetry) 2522 2523 IF (almo_scf_env%deloc_truncate_virt .NE. virt_full) THEN 2524 CALL matrix_almo_create(matrix_new=almo_scf_env%opt_k_t_rr(ispin), & 2525 matrix_qs=matrix_s0, & 2526 almo_scf_env=almo_scf_env, & 2527 name_new="OPT_K_U_RR", & 2528 size_keys=(/almo_mat_dim_virt, almo_mat_dim_virt/), & 2529 symmetry_new=dbcsr_type_no_symmetry, & 2530 spin_key=ispin, & 2531 init_domains=.FALSE.) 2532 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_vv_disc(ispin), & 2533 matrix_qs=matrix_s0, & 2534 almo_scf_env=almo_scf_env, & 2535 name_new="VV_DISC", & 2536 size_keys=(/almo_mat_dim_virt_disc, almo_mat_dim_virt_disc/), & 2537 symmetry_new=dbcsr_type_symmetric, & 2538 spin_key=ispin, & 2539 init_domains=.FALSE.) 2540 CALL matrix_almo_create(matrix_new=almo_scf_env%opt_k_t_dd(ispin), & 2541 matrix_qs=matrix_s0, & 2542 almo_scf_env=almo_scf_env, & 2543 name_new="OPT_K_U_DD", & 2544 size_keys=(/almo_mat_dim_virt_disc, almo_mat_dim_virt_disc/), & 2545 symmetry_new=dbcsr_type_no_symmetry, & 2546 spin_key=ispin, & 2547 init_domains=.FALSE.) 2548 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_vv_disc_blk(ispin), & 2549 matrix_qs=matrix_s0, & 2550 almo_scf_env=almo_scf_env, & 2551 name_new="VV_DISC_BLK", & 2552 size_keys=(/almo_mat_dim_virt_disc, almo_mat_dim_virt_disc/), & 2553 symmetry_new=dbcsr_type_symmetric, & 2554 spin_key=ispin, & 2555 init_domains=.TRUE.) 2556 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_k_blk(ispin), & 2557 matrix_qs=matrix_s0, & 2558 almo_scf_env=almo_scf_env, & 2559 name_new="K_BLK", & 2560 size_keys=(/almo_mat_dim_virt_disc, almo_mat_dim_virt/), & 2561 symmetry_new=dbcsr_type_no_symmetry, & 2562 spin_key=ispin, & 2563 init_domains=.TRUE.) 2564 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_k_blk_ones(ispin), & 2565 matrix_qs=matrix_s0, & 2566 almo_scf_env=almo_scf_env, & 2567 name_new="K_BLK_1", & 2568 size_keys=(/almo_mat_dim_virt_disc, almo_mat_dim_virt/), & 2569 symmetry_new=dbcsr_type_no_symmetry, & 2570 spin_key=ispin, & 2571 init_domains=.TRUE.) 2572 CALL matrix_almo_create(matrix_new=almo_scf_env%opt_k_denom(ispin), & 2573 matrix_qs=matrix_s0, & 2574 almo_scf_env=almo_scf_env, & 2575 name_new="OPT_K_DENOM", & 2576 size_keys=(/almo_mat_dim_virt_disc, almo_mat_dim_virt/), & 2577 symmetry_new=dbcsr_type_no_symmetry, & 2578 spin_key=ispin, & 2579 init_domains=.FALSE.) 2580 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_k_tr(ispin), & 2581 matrix_qs=matrix_s0, & 2582 almo_scf_env=almo_scf_env, & 2583 name_new="K_TR", & 2584 size_keys=(/almo_mat_dim_virt, almo_mat_dim_virt_disc/), & 2585 symmetry_new=dbcsr_type_no_symmetry, & 2586 spin_key=ispin, & 2587 init_domains=.FALSE.) 2588 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_v_disc_blk(ispin), & 2589 matrix_qs=matrix_s0, & 2590 almo_scf_env=almo_scf_env, & 2591 name_new="V_DISC_BLK", & 2592 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_virt_disc/), & 2593 symmetry_new=dbcsr_type_no_symmetry, & 2594 spin_key=ispin, & 2595 init_domains=.FALSE.) 2596 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_v_disc(ispin), & 2597 matrix_qs=matrix_s0, & 2598 almo_scf_env=almo_scf_env, & 2599 name_new="V_DISC", & 2600 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_virt_disc/), & 2601 symmetry_new=dbcsr_type_no_symmetry, & 2602 spin_key=ispin, & 2603 init_domains=.FALSE.) 2604 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_ov_disc(ispin), & 2605 matrix_qs=matrix_s0, & 2606 almo_scf_env=almo_scf_env, & 2607 name_new="OV_DISC", & 2608 size_keys=(/almo_mat_dim_occ, almo_mat_dim_virt_disc/), & 2609 symmetry_new=dbcsr_type_no_symmetry, & 2610 spin_key=ispin, & 2611 init_domains=.FALSE.) 2612 2613 ENDIF ! end need_discarded_virtuals 2614 2615 ENDDO ! spin 2616 ENDIF 2617 2618 ! create matrices of orbital energies if necessary 2619 IF (almo_scf_env%need_orbital_energies) THEN 2620 ALLOCATE (almo_scf_env%matrix_eoo(nspins)) 2621 ALLOCATE (almo_scf_env%matrix_evv_full(nspins)) 2622 DO ispin = 1, nspins 2623 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_eoo(ispin), & 2624 matrix_qs=matrix_s0, & 2625 almo_scf_env=almo_scf_env, & 2626 name_new="E_OCC", & 2627 size_keys=(/almo_mat_dim_occ, almo_mat_dim_occ/), & 2628 symmetry_new=dbcsr_type_no_symmetry, & 2629 spin_key=ispin, & 2630 init_domains=.FALSE.) 2631 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_evv_full(ispin), & 2632 matrix_qs=matrix_s0, & 2633 almo_scf_env=almo_scf_env, & 2634 name_new="E_VIRT", & 2635 size_keys=(/almo_mat_dim_virt_full, almo_mat_dim_virt_full/), & 2636 symmetry_new=dbcsr_type_no_symmetry, & 2637 spin_key=ispin, & 2638 init_domains=.FALSE.) 2639 ENDDO 2640 ENDIF 2641 2642 ! Density and KS matrices 2643 ALLOCATE (almo_scf_env%matrix_p(nspins)) 2644 ALLOCATE (almo_scf_env%matrix_p_blk(nspins)) 2645 ALLOCATE (almo_scf_env%matrix_ks(nspins)) 2646 ALLOCATE (almo_scf_env%matrix_ks_blk(nspins)) 2647 IF (almo_scf_env%need_previous_ks) & 2648 ALLOCATE (almo_scf_env%matrix_ks_0deloc(nspins)) 2649 DO ispin = 1, nspins 2650 ! RZK-warning copy with symmery but remember that this might cause problems 2651 CALL dbcsr_create(almo_scf_env%matrix_p(ispin), & 2652 template=almo_scf_env%matrix_s(1), & 2653 matrix_type=dbcsr_type_symmetric) 2654 CALL dbcsr_create(almo_scf_env%matrix_ks(ispin), & 2655 template=almo_scf_env%matrix_s(1), & 2656 matrix_type=dbcsr_type_symmetric) 2657 IF (almo_scf_env%need_previous_ks) THEN 2658 CALL dbcsr_create(almo_scf_env%matrix_ks_0deloc(ispin), & 2659 template=almo_scf_env%matrix_s(1), & 2660 matrix_type=dbcsr_type_symmetric) 2661 ENDIF 2662 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_p_blk(ispin), & 2663 matrix_qs=matrix_s0, & 2664 almo_scf_env=almo_scf_env, & 2665 name_new="P_BLK", & 2666 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_aobasis/), & 2667 symmetry_new=dbcsr_type_symmetric, & 2668 spin_key=ispin, & 2669 init_domains=.TRUE.) 2670 CALL matrix_almo_create(matrix_new=almo_scf_env%matrix_ks_blk(ispin), & 2671 matrix_qs=matrix_s0, & 2672 almo_scf_env=almo_scf_env, & 2673 name_new="KS_BLK", & 2674 size_keys=(/almo_mat_dim_aobasis, almo_mat_dim_aobasis/), & 2675 symmetry_new=dbcsr_type_symmetric, & 2676 spin_key=ispin, & 2677 init_domains=.TRUE.) 2678 ENDDO 2679 2680 CALL timestop(handle) 2681 2682 END SUBROUTINE almo_scf_env_create_matrices 2683 2684! ************************************************************************************************** 2685!> \brief clean up procedures for almo scf 2686!> \param almo_scf_env ... 2687!> \par History 2688!> 2011.06 created [Rustam Z Khaliullin] 2689!> 2018.09 smearing support [Ruben Staub] 2690!> \author Rustam Z Khaliullin 2691! ************************************************************************************************** 2692 SUBROUTINE almo_scf_clean_up(almo_scf_env) 2693 2694 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env 2695 2696 CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_clean_up', & 2697 routineP = moduleN//':'//routineN 2698 2699 INTEGER :: handle, ispin, unit_nr 2700 TYPE(cp_logger_type), POINTER :: logger 2701 2702 CALL timeset(routineN, handle) 2703 2704 ! get a useful output_unit 2705 logger => cp_get_default_logger() 2706 IF (logger%para_env%ionode) THEN 2707 unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.) 2708 ELSE 2709 unit_nr = -1 2710 ENDIF 2711 2712 ! release matrices 2713 CALL dbcsr_release(almo_scf_env%matrix_s(1)) 2714 CALL dbcsr_release(almo_scf_env%matrix_s_blk(1)) 2715 IF (almo_scf_env%almo_update_algorithm .EQ. almo_scf_diag) THEN 2716 CALL dbcsr_release(almo_scf_env%matrix_s_blk_sqrt_inv(1)) 2717 CALL dbcsr_release(almo_scf_env%matrix_s_blk_sqrt(1)) 2718 ELSE IF (almo_scf_env%almo_update_algorithm .EQ. almo_scf_dm_sign) THEN 2719 CALL dbcsr_release(almo_scf_env%matrix_s_blk_inv(1)) 2720 ENDIF 2721 DO ispin = 1, almo_scf_env%nspins 2722 CALL dbcsr_release(almo_scf_env%quench_t(ispin)) 2723 CALL dbcsr_release(almo_scf_env%quench_t_blk(ispin)) 2724 CALL dbcsr_release(almo_scf_env%matrix_t_blk(ispin)) 2725 CALL dbcsr_release(almo_scf_env%matrix_err_blk(ispin)) 2726 CALL dbcsr_release(almo_scf_env%matrix_err_xx(ispin)) 2727 CALL dbcsr_release(almo_scf_env%matrix_t_tr(ispin)) 2728 CALL dbcsr_release(almo_scf_env%matrix_sigma(ispin)) 2729 CALL dbcsr_release(almo_scf_env%matrix_sigma_blk(ispin)) 2730 CALL dbcsr_release(almo_scf_env%matrix_sigma_inv_0deloc(ispin)) 2731 CALL dbcsr_release(almo_scf_env%matrix_sigma_inv(ispin)) 2732 CALL dbcsr_release(almo_scf_env%matrix_t(ispin)) 2733 CALL dbcsr_release(almo_scf_env%matrix_sigma_sqrt(ispin)) 2734 CALL dbcsr_release(almo_scf_env%matrix_sigma_sqrt_inv(ispin)) 2735 CALL dbcsr_release(almo_scf_env%matrix_p(ispin)) 2736 CALL dbcsr_release(almo_scf_env%matrix_ks(ispin)) 2737 CALL dbcsr_release(almo_scf_env%matrix_p_blk(ispin)) 2738 CALL dbcsr_release(almo_scf_env%matrix_ks_blk(ispin)) 2739 IF (almo_scf_env%need_previous_ks) THEN 2740 CALL dbcsr_release(almo_scf_env%matrix_ks_0deloc(ispin)) 2741 ENDIF 2742 IF (almo_scf_env%need_virtuals) THEN 2743 CALL dbcsr_release(almo_scf_env%matrix_v_blk(ispin)) 2744 CALL dbcsr_release(almo_scf_env%matrix_v_full_blk(ispin)) 2745 CALL dbcsr_release(almo_scf_env%matrix_v(ispin)) 2746 CALL dbcsr_release(almo_scf_env%matrix_vo(ispin)) 2747 CALL dbcsr_release(almo_scf_env%matrix_x(ispin)) 2748 CALL dbcsr_release(almo_scf_env%matrix_ov(ispin)) 2749 CALL dbcsr_release(almo_scf_env%matrix_ov_full(ispin)) 2750 CALL dbcsr_release(almo_scf_env%matrix_sigma_vv(ispin)) 2751 CALL dbcsr_release(almo_scf_env%matrix_sigma_vv_blk(ispin)) 2752 CALL dbcsr_release(almo_scf_env%matrix_sigma_vv_sqrt(ispin)) 2753 CALL dbcsr_release(almo_scf_env%matrix_sigma_vv_sqrt_inv(ispin)) 2754 CALL dbcsr_release(almo_scf_env%matrix_vv_full_blk(ispin)) 2755 IF (almo_scf_env%deloc_truncate_virt .NE. virt_full) THEN 2756 CALL dbcsr_release(almo_scf_env%matrix_k_tr(ispin)) 2757 CALL dbcsr_release(almo_scf_env%matrix_k_blk(ispin)) 2758 CALL dbcsr_release(almo_scf_env%matrix_k_blk_ones(ispin)) 2759 CALL dbcsr_release(almo_scf_env%matrix_v_disc(ispin)) 2760 CALL dbcsr_release(almo_scf_env%matrix_v_disc_blk(ispin)) 2761 CALL dbcsr_release(almo_scf_env%matrix_ov_disc(ispin)) 2762 CALL dbcsr_release(almo_scf_env%matrix_vv_disc_blk(ispin)) 2763 CALL dbcsr_release(almo_scf_env%matrix_vv_disc(ispin)) 2764 CALL dbcsr_release(almo_scf_env%opt_k_t_dd(ispin)) 2765 CALL dbcsr_release(almo_scf_env%opt_k_t_rr(ispin)) 2766 CALL dbcsr_release(almo_scf_env%opt_k_denom(ispin)) 2767 ENDIF 2768 ENDIF 2769 IF (almo_scf_env%need_orbital_energies) THEN 2770 CALL dbcsr_release(almo_scf_env%matrix_eoo(ispin)) 2771 CALL dbcsr_release(almo_scf_env%matrix_evv_full(ispin)) 2772 ENDIF 2773 ENDDO 2774 2775 ! deallocate matrices 2776 DEALLOCATE (almo_scf_env%matrix_p) 2777 DEALLOCATE (almo_scf_env%matrix_p_blk) 2778 DEALLOCATE (almo_scf_env%matrix_ks) 2779 DEALLOCATE (almo_scf_env%matrix_ks_blk) 2780 DEALLOCATE (almo_scf_env%matrix_t_blk) 2781 DEALLOCATE (almo_scf_env%matrix_err_blk) 2782 DEALLOCATE (almo_scf_env%matrix_err_xx) 2783 DEALLOCATE (almo_scf_env%matrix_t) 2784 DEALLOCATE (almo_scf_env%matrix_t_tr) 2785 DEALLOCATE (almo_scf_env%matrix_sigma) 2786 DEALLOCATE (almo_scf_env%matrix_sigma_blk) 2787 DEALLOCATE (almo_scf_env%matrix_sigma_inv_0deloc) 2788 DEALLOCATE (almo_scf_env%matrix_sigma_sqrt) 2789 DEALLOCATE (almo_scf_env%matrix_sigma_sqrt_inv) 2790 DEALLOCATE (almo_scf_env%matrix_sigma_inv) 2791 DEALLOCATE (almo_scf_env%quench_t) 2792 DEALLOCATE (almo_scf_env%quench_t_blk) 2793 IF (almo_scf_env%need_virtuals) THEN 2794 DEALLOCATE (almo_scf_env%matrix_v_blk) 2795 DEALLOCATE (almo_scf_env%matrix_v_full_blk) 2796 DEALLOCATE (almo_scf_env%matrix_v) 2797 DEALLOCATE (almo_scf_env%matrix_vo) 2798 DEALLOCATE (almo_scf_env%matrix_x) 2799 DEALLOCATE (almo_scf_env%matrix_ov) 2800 DEALLOCATE (almo_scf_env%matrix_ov_full) 2801 DEALLOCATE (almo_scf_env%matrix_sigma_vv) 2802 DEALLOCATE (almo_scf_env%matrix_sigma_vv_blk) 2803 DEALLOCATE (almo_scf_env%matrix_sigma_vv_sqrt) 2804 DEALLOCATE (almo_scf_env%matrix_sigma_vv_sqrt_inv) 2805 DEALLOCATE (almo_scf_env%matrix_vv_full_blk) 2806 IF (almo_scf_env%deloc_truncate_virt .NE. virt_full) THEN 2807 DEALLOCATE (almo_scf_env%matrix_k_tr) 2808 DEALLOCATE (almo_scf_env%matrix_k_blk) 2809 DEALLOCATE (almo_scf_env%matrix_v_disc) 2810 DEALLOCATE (almo_scf_env%matrix_v_disc_blk) 2811 DEALLOCATE (almo_scf_env%matrix_ov_disc) 2812 DEALLOCATE (almo_scf_env%matrix_vv_disc_blk) 2813 DEALLOCATE (almo_scf_env%matrix_vv_disc) 2814 DEALLOCATE (almo_scf_env%matrix_k_blk_ones) 2815 DEALLOCATE (almo_scf_env%opt_k_t_dd) 2816 DEALLOCATE (almo_scf_env%opt_k_t_rr) 2817 DEALLOCATE (almo_scf_env%opt_k_denom) 2818 ENDIF 2819 ENDIF 2820 IF (almo_scf_env%need_previous_ks) THEN 2821 DEALLOCATE (almo_scf_env%matrix_ks_0deloc) 2822 ENDIF 2823 IF (almo_scf_env%need_orbital_energies) THEN 2824 DEALLOCATE (almo_scf_env%matrix_eoo) 2825 DEALLOCATE (almo_scf_env%matrix_evv_full) 2826 ENDIF 2827 2828 ! clean up other variables 2829 DO ispin = 1, almo_scf_env%nspins 2830 CALL release_submatrices( & 2831 almo_scf_env%domain_preconditioner(:, ispin)) 2832 CALL release_submatrices(almo_scf_env%domain_s_inv(:, ispin)) 2833 CALL release_submatrices(almo_scf_env%domain_s_sqrt_inv(:, ispin)) 2834 CALL release_submatrices(almo_scf_env%domain_s_sqrt(:, ispin)) 2835 CALL release_submatrices(almo_scf_env%domain_ks_xx(:, ispin)) 2836 CALL release_submatrices(almo_scf_env%domain_t(:, ispin)) 2837 CALL release_submatrices(almo_scf_env%domain_err(:, ispin)) 2838 CALL release_submatrices(almo_scf_env%domain_r_down_up(:, ispin)) 2839 ENDDO 2840 DEALLOCATE (almo_scf_env%domain_preconditioner) 2841 DEALLOCATE (almo_scf_env%domain_s_inv) 2842 DEALLOCATE (almo_scf_env%domain_s_sqrt_inv) 2843 DEALLOCATE (almo_scf_env%domain_s_sqrt) 2844 DEALLOCATE (almo_scf_env%domain_ks_xx) 2845 DEALLOCATE (almo_scf_env%domain_t) 2846 DEALLOCATE (almo_scf_env%domain_err) 2847 DEALLOCATE (almo_scf_env%domain_r_down_up) 2848 DO ispin = 1, almo_scf_env%nspins 2849 DEALLOCATE (almo_scf_env%domain_map(ispin)%pairs) 2850 DEALLOCATE (almo_scf_env%domain_map(ispin)%index1) 2851 ENDDO 2852 DEALLOCATE (almo_scf_env%domain_map) 2853 DEALLOCATE (almo_scf_env%domain_index_of_ao) 2854 DEALLOCATE (almo_scf_env%domain_index_of_atom) 2855 DEALLOCATE (almo_scf_env%first_atom_of_domain) 2856 DEALLOCATE (almo_scf_env%last_atom_of_domain) 2857 DEALLOCATE (almo_scf_env%nbasis_of_domain) 2858 DEALLOCATE (almo_scf_env%nocc_of_domain) 2859 DEALLOCATE (almo_scf_env%real_ne_of_domain) 2860 DEALLOCATE (almo_scf_env%nvirt_full_of_domain) 2861 DEALLOCATE (almo_scf_env%nvirt_of_domain) 2862 DEALLOCATE (almo_scf_env%nvirt_disc_of_domain) 2863 DEALLOCATE (almo_scf_env%mu_of_domain) 2864 DEALLOCATE (almo_scf_env%cpu_of_domain) 2865 DEALLOCATE (almo_scf_env%charge_of_domain) 2866 DEALLOCATE (almo_scf_env%multiplicity_of_domain) 2867 IF (almo_scf_env%smear) THEN 2868 DEALLOCATE (almo_scf_env%mo_energies) 2869 DEALLOCATE (almo_scf_env%kTS) 2870 END IF 2871 2872 DEALLOCATE (almo_scf_env%domain_index_of_ao_block) 2873 DEALLOCATE (almo_scf_env%domain_index_of_mo_block) 2874 2875 CALL cp_para_env_release(almo_scf_env%para_env) 2876 CALL cp_blacs_env_release(almo_scf_env%blacs_env) 2877 2878 CALL timestop(handle) 2879 2880 END SUBROUTINE almo_scf_clean_up 2881 2882! ************************************************************************************************** 2883!> \brief Do post scf calculations with ALMO 2884!> WARNING: ALMO post scf calculation may not work for certain quantities, 2885!> like forces, since ALMO wave function is only 'partially' optimized 2886!> \param qs_env ... 2887!> \par History 2888!> 2016.12 created [Yifei Shi] 2889!> \author Yifei Shi 2890! ************************************************************************************************** 2891 SUBROUTINE almo_post_scf_compute_properties(qs_env) 2892 TYPE(qs_environment_type), POINTER :: qs_env 2893 2894 CHARACTER(len=*), PARAMETER :: routineN = 'almo_post_scf_compute_properties', & 2895 routineP = moduleN//':'//routineN 2896 2897 CALL qs_scf_compute_properties(qs_env) 2898 2899 END SUBROUTINE almo_post_scf_compute_properties 2900 2901END MODULE almo_scf 2902 2903