comp-chem/cp2k
1
0
Fork 0
mirror of https://github.com/cp2k/cp2k.git synced 2026-07-28 22:25:32 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

ALMO Update (Rustam Khalliulin)

Browse source 
* ALMO SCF for strongly-interacting fragments
* Langevin AIMD based on under-converged ALMOs for weakly-interacting fragments
* ALMO code for zero-electron fragments and charged fragments
* wider variety of ALMO options
* major refactoring of the ALMO code for the current and future developments
* some new regtests, some changes in regtest values


svn-origin-rev: 18489
This commit is contained in:
Jürg Hutter 2018-06-26 15:22:04 +00:00
parent e0901beed2
commit b28bda4cfa
36 changed files with 6158 additions and 2782 deletions
Download patch file Download diff file Expand all files Collapse all files

4
src/admm_dm_methods.F
View file

@ -135,6 +135,7 @@ CONTAINS
INTEGER :: ispin
LOGICAL :: s_mstruct_changed
REAL(KIND=dp) :: threshold
TYPE(admm_dm_type), POINTER :: admm_dm
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s_aux, matrix_s_mixed, rho_ao, &
rho_ao_aux
@ -160,7 +161,8 @@ CONTAINS
IF (s_mstruct_changed) THEN
! Calculate A = S_aux^(-1) * S_mixed
CALL dbcsr_create(matrix_s_aux_inv, template=matrix_s_aux(1)%matrix, matrix_type="N")
CALL invert_Hotelling(matrix_s_aux_inv, matrix_s_aux(1)%matrix, admm_dm%eps_filter)
threshold = MAX(admm_dm%eps_filter, 1.0e-12_dp)
CALL invert_Hotelling(matrix_s_aux_inv, matrix_s_aux(1)%matrix, threshold)
IF (.NOT. ASSOCIATED(admm_dm%matrix_A)) THEN
ALLOCATE (admm_dm%matrix_A)

407
src/almo_scf.F
View file

@ -12,18 +12,18 @@
! **************************************************************************************************
MODULE almo_scf
USE almo_scf_methods, ONLY: almo_scf_p_blk_to_t_blk,&
almo_scf_t_blk_to_p,&
almo_scf_t_blk_to_t_blk_orthonormal,&
almo_scf_t_to_proj,&
distribute_domains,&
orthogonalize_mos
USE almo_scf_optimizer, ONLY: almo_scf_block_diagonal,&
almo_scf_xalmo_eigensolver,&
almo_scf_xalmo_pcg
USE almo_scf_qs, ONLY: almo_scf_construct_quencher,&
almo_scf_init_qs,&
USE almo_scf_qs, ONLY: almo_dm_to_almo_ks,&
almo_scf_construct_quencher,&
calculate_w_matrix_almo,&
construct_qs_mos,&
init_almo_ks_matrix_via_qs,&
matrix_almo_create,&
matrix_almo_to_qs,&
matrix_qs_to_almo
USE almo_scf_types, ONLY: almo_mat_dim_aobasis,&
almo_mat_dim_occ,&
@ -36,6 +36,7 @@ MODULE almo_scf
USE atomic_kind_types, ONLY: atomic_kind_type
USE bibliography, ONLY: Khaliullin2013,&
Kolafa2004,&
Scheiber2018,&
cite_reference
USE cp_blacs_env, ONLY: cp_blacs_env_release,&
cp_blacs_env_retain
@ -60,9 +61,10 @@ MODULE almo_scf
almo_deloc_none, almo_deloc_scf, almo_deloc_x, almo_deloc_x_then_scf, &
almo_deloc_xalmo_1diag, almo_deloc_xalmo_scf, almo_deloc_xalmo_x, almo_deloc_xk, &
almo_domain_layout_molecular, almo_mat_distr_atomic, almo_mat_distr_molecular, &
almo_scf_diag, almo_scf_dm_sign, almo_scf_pcg, atomic_guess, molecular_guess, &
optimizer_diis, optimizer_pcg, restart_guess, virt_full, virt_number, virt_occ_size, &
xalmo_case_block_diag, xalmo_case_fully_deloc, xalmo_case_normal
almo_scf_diag, almo_scf_dm_sign, almo_scf_pcg, almo_scf_skip, atomic_guess, &
molecular_guess, optimizer_diis, optimizer_lin_eq_pcg, optimizer_pcg, restart_guess, &
virt_full, virt_number, virt_occ_size, xalmo_case_block_diag, xalmo_case_fully_deloc, &
xalmo_case_normal, xalmo_trial_r0_out
USE input_section_types, ONLY: section_vals_get_subs_vals,&
section_vals_type
USE iterate_matrix, ONLY: invert_Hotelling,&
@ -75,17 +77,13 @@ MODULE almo_scf
USE mscfg_types, ONLY: get_matrix_from_submatrices,&
molecular_scf_guess_env_type
USE particle_types, ONLY: particle_type
USE qs_energy_types, ONLY: qs_energy_type
USE qs_environment_types, ONLY: get_qs_env,&
qs_environment_type
USE qs_initial_guess, ONLY: calculate_atomic_block_dm,&
calculate_mopac_dm
USE qs_kind_types, ONLY: qs_kind_type
USE qs_ks_methods, ONLY: qs_ks_update_qs_env
USE qs_ks_types, ONLY: qs_ks_did_change
USE qs_mo_types, ONLY: get_mo_set,&
mo_set_p_type
USE qs_rho_methods, ONLY: qs_rho_update_rho
USE qs_rho_types, ONLY: qs_rho_get,&
qs_rho_type
USE qs_scf_post_scf, ONLY: qs_scf_compute_properties
@ -196,6 +194,7 @@ CONTAINS
almo_scf_env%opt_block_diag_pcg%optimizer_type = optimizer_pcg
almo_scf_env%opt_xalmo_diis%optimizer_type = optimizer_diis
almo_scf_env%opt_xalmo_pcg%optimizer_type = optimizer_pcg
almo_scf_env%opt_xalmo_newton_pcg_solver%optimizer_type = optimizer_lin_eq_pcg
! get info from the qs_env
CALL get_qs_env(qs_env, &
@ -223,9 +222,6 @@ CONTAINS
nmols = almo_scf_env%nmolecules
natoms = almo_scf_env%natoms
! parse the almo_scf section and set appropriate quantities
!CALL almo_scf_init_read_write_input(input, almo_scf_env)
! Define groups: either atomic or molecular
IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular) THEN
almo_scf_env%ndomains = almo_scf_env%nmolecules
@ -366,22 +362,57 @@ CONTAINS
almo_scf_env%calc_forces = calc_forces
IF (calc_forces) THEN
IF (almo_scf_env%deloc_method == almo_deloc_x .OR. &
almo_scf_env%deloc_method == almo_deloc_xalmo_x .OR. &
almo_scf_env%deloc_method == almo_deloc_xalmo_1diag) &
CPABORT(" Forces are not implemented for this method yet ")
ENDIF
almo_scf_env%calc_forces = calc_forces
IF (calc_forces) THEN
IF (almo_scf_env%deloc_method == almo_deloc_x .OR. &
almo_scf_env%deloc_method == almo_deloc_xalmo_x .OR. &
almo_scf_env%deloc_method == almo_deloc_xalmo_1diag) &
CPABORT(" Forces are not implemented for this method yet ")
CALL cite_reference(Scheiber2018)
IF (almo_scf_env%deloc_method .EQ. almo_deloc_x .OR. &
almo_scf_env%deloc_method .EQ. almo_deloc_xalmo_x .OR. &
almo_scf_env%deloc_method .EQ. almo_deloc_xalmo_1diag) THEN
CPABORT("Forces for perturbative methods are NYI. Change DELOCALIZE_METHOD")
ENDIF
! switch to ASPC after a certain number of exact steps is done
IF (almo_scf_env%almo_history%istore .GT. (almo_scf_env%almo_history%nstore+1)) THEN
IF (almo_scf_env%opt_block_diag_pcg%eps_error_early .GT. 0.0_dp) THEN
almo_scf_env%opt_block_diag_pcg%eps_error = almo_scf_env%opt_block_diag_pcg%eps_error_early
almo_scf_env%opt_block_diag_pcg%early_stopping_on = .TRUE.
IF (unit_nr > 0) WRITE (*, *) "ALMO_OPTIMIZER_PCG: EPS_ERROR_EARLY is on"
ENDIF
IF (almo_scf_env%opt_block_diag_diis%eps_error_early .GT. 0.0_dp) THEN
almo_scf_env%opt_block_diag_diis%eps_error = almo_scf_env%opt_block_diag_diis%eps_error_early
almo_scf_env%opt_block_diag_diis%early_stopping_on = .TRUE.
IF (unit_nr > 0) WRITE (*, *) "ALMO_OPTIMIZER_DIIS: EPS_ERROR_EARLY is on"
ENDIF
IF (almo_scf_env%opt_block_diag_pcg%max_iter_early .GT. 0) THEN
almo_scf_env%opt_block_diag_pcg%max_iter = almo_scf_env%opt_block_diag_pcg%max_iter_early
almo_scf_env%opt_block_diag_pcg%early_stopping_on = .TRUE.
IF (unit_nr > 0) WRITE (*, *) "ALMO_OPTIMIZER_PCG: MAX_ITER_EARLY is on"
ENDIF
IF (almo_scf_env%opt_block_diag_diis%max_iter_early .GT. 0) THEN
almo_scf_env%opt_block_diag_diis%max_iter = almo_scf_env%opt_block_diag_diis%max_iter_early
almo_scf_env%opt_block_diag_diis%early_stopping_on = .TRUE.
IF (unit_nr > 0) WRITE (*, *) "ALMO_OPTIMIZER_DIIS: MAX_ITER_EARLY is on"
ENDIF
ELSE
almo_scf_env%opt_block_diag_diis%early_stopping_on = .FALSE.
almo_scf_env%opt_block_diag_pcg%early_stopping_on = .FALSE.
ENDIF
IF (almo_scf_env%xalmo_history%istore .GT. (almo_scf_env%xalmo_history%nstore+1)) THEN
IF (almo_scf_env%opt_xalmo_pcg%eps_error_early .GT. 0.0_dp) THEN
almo_scf_env%opt_xalmo_pcg%eps_error = almo_scf_env%opt_xalmo_pcg%eps_error_early
almo_scf_env%opt_xalmo_pcg%early_stopping_on = .TRUE.
IF (unit_nr > 0) WRITE (*, *) "XALMO_OPTIMIZER_PCG: EPS_ERROR_EARLY is on"
ENDIF
IF (almo_scf_env%opt_xalmo_pcg%max_iter_early .GT. 0.0_dp) THEN
almo_scf_env%opt_xalmo_pcg%max_iter = almo_scf_env%opt_xalmo_pcg%max_iter_early
almo_scf_env%opt_xalmo_pcg%early_stopping_on = .TRUE.
IF (unit_nr > 0) WRITE (*, *) "XALMO_OPTIMIZER_PCG: MAX_ITER_EARLY is on"
ENDIF
ELSE
almo_scf_env%opt_xalmo_pcg%early_stopping_on = .FALSE.
ENDIF
ENDIF
! create all matrices
CALL almo_scf_env_create_matrices(almo_scf_env, matrix_s(1)%matrix)
! set up matrix S and all required functions of S
almo_scf_env%s_inv_done = .FALSE.
almo_scf_env%s_sqrt_done = .FALSE.
@ -416,8 +447,12 @@ CONTAINS
ALLOCATE (almo_scf_env%domain_r_down_up(ndomains, nspins))
CALL init_submatrices(almo_scf_env%domain_r_down_up)
! initialization of the QS settings with the ALMO flavor
CALL almo_scf_init_qs(qs_env, almo_scf_env)
! initialization of the KS matrix
CALL init_almo_ks_matrix_via_qs(qs_env, &
almo_scf_env%matrix_ks, &
almo_scf_env%mat_distr_aos, &
almo_scf_env%eps_filter)
CALL construct_qs_mos(qs_env, almo_scf_env)
CALL timestop(handle)
@ -443,7 +478,7 @@ CONTAINS
nspins, unit_nr
INTEGER, DIMENSION(2) :: nelectron_spin
LOGICAL :: aspc_guess, has_unit_metric
REAL(KIND=dp) :: alpha, cs_pos
REAL(KIND=dp) :: alpha, cs_pos, energy
TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
TYPE(cp_logger_type), POINTER :: logger
TYPE(cp_para_env_type), POINTER :: para_env
@ -452,7 +487,6 @@ CONTAINS
TYPE(dft_control_type), POINTER :: dft_control
TYPE(molecular_scf_guess_env_type), POINTER :: mscfg_env
TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
TYPE(qs_energy_type), POINTER :: energy
TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
TYPE(qs_rho_type), POINTER :: rho
@ -472,7 +506,6 @@ CONTAINS
CALL get_qs_env(qs_env, &
dft_control=dft_control, &
matrix_s=matrix_s, &
energy=energy, &
atomic_kind_set=atomic_kind_set, &
qs_kind_set=qs_kind_set, &
particle_set=particle_set, &
@ -532,7 +565,7 @@ CONTAINS
DO ispin = 1, nspins
! copy the atomic-block dm into matrix_p_blk
CALL matrix_qs_to_almo(rho_ao(ispin)%matrix, &
almo_scf_env%matrix_p_blk(ispin), almo_scf_env, &
almo_scf_env%matrix_p_blk(ispin), almo_scf_env%mat_distr_aos, &
.FALSE.)
CALL dbcsr_filter(almo_scf_env%matrix_p_blk(ispin), &
almo_scf_env%eps_filter)
@ -597,27 +630,55 @@ CONTAINS
ENDIF !aspc_guess?
CALL almo_scf_t_blk_to_t_blk_orthonormal(almo_scf_env)
CALL almo_scf_t_blk_to_p(almo_scf_env, &
use_sigma_inv_guess=.FALSE.)
DO ispin = 1, nspins
CALL matrix_almo_to_qs(almo_scf_env%matrix_p(ispin), &
rho_ao(ispin)%matrix, &
almo_scf_env)
CALL orthogonalize_mos(ket=almo_scf_env%matrix_t_blk(ispin), &
overlap=almo_scf_env%matrix_sigma_blk(ispin), &
metric=almo_scf_env%matrix_s_blk(1), &
retain_locality=.TRUE., &
only_normalize=.FALSE., &
nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), &
eps_filter=almo_scf_env%eps_filter, &
order_lanczos=almo_scf_env%order_lanczos, &
eps_lanczos=almo_scf_env%eps_lanczos, &
max_iter_lanczos=almo_scf_env%max_iter_lanczos)
CALL almo_scf_t_to_proj(t=almo_scf_env%matrix_t_blk(ispin), &
p=almo_scf_env%matrix_p(ispin), &
eps_filter=almo_scf_env%eps_filter, &
orthog_orbs=.FALSE., &
nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), &
s=almo_scf_env%matrix_s(1), &
sigma=almo_scf_env%matrix_sigma(ispin), &
sigma_inv=almo_scf_env%matrix_sigma_inv(ispin), &
use_guess=.FALSE., &
algorithm=almo_scf_env%sigma_inv_algorithm, &
eps_lanczos=almo_scf_env%eps_lanczos, &
max_iter_lanczos=almo_scf_env%max_iter_lanczos, &
inv_eps_factor=almo_scf_env%matrix_iter_eps_error_factor, &
para_env=almo_scf_env%para_env, &
blacs_env=almo_scf_env%blacs_env)
ENDDO
CALL qs_rho_update_rho(rho, qs_env=qs_env)
CALL qs_ks_did_change(qs_env%ks_env, rho_changed=.TRUE.)
CALL qs_ks_update_qs_env(qs_env, calculate_forces=.FALSE., &
just_energy=.FALSE.)
! compute dm from the projector(s)
IF (nspins == 1) THEN
CALL dbcsr_scale(almo_scf_env%matrix_p(1), 2.0_dp)
ENDIF
CALL almo_dm_to_almo_ks(qs_env, &
almo_scf_env%matrix_p, &
almo_scf_env%matrix_ks, &
energy, &
almo_scf_env%eps_filter, &
almo_scf_env%mat_distr_aos)
IF (unit_nr > 0) THEN
IF (almo_scf_env%almo_scf_guess .EQ. molecular_guess) THEN
WRITE (unit_nr, '(T2,A38,F40.10)') "Single-molecule energy:", &
SUM(mscfg_env%energy_of_frag)
ENDIF
WRITE (unit_nr, '(T2,A38,F40.10)') "Energy of the initial guess:", energy%total
WRITE (unit_nr, '(T2,A38,F40.10)') "Energy of the initial guess:", energy
WRITE (unit_nr, '()')
ENDIF
@ -632,10 +693,10 @@ CONTAINS
!> 2016.11 created [Rustam Z Khaliullin]
!> \author Rustam Khaliullin
! **************************************************************************************************
SUBROUTINE almo_scf_store_result(almo_scf_env)
SUBROUTINE almo_scf_store_extrapolation_data(almo_scf_env)
TYPE(almo_scf_env_type) :: almo_scf_env
CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_store_result', &
CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_store_extrapolation_data', &
routineP = moduleN//':'//routineN
INTEGER :: handle, ispin, istore, unit_nr
@ -769,7 +830,7 @@ CONTAINS
CALL timestop(handle)
END SUBROUTINE almo_scf_store_result
END SUBROUTINE almo_scf_store_extrapolation_data
! **************************************************************************************************
!> \brief Prints out a short summary about the ALMO SCF job
@ -789,7 +850,8 @@ CONTAINS
CHARACTER(len=13) :: neig_string
CHARACTER(len=33) :: deloc_method_string
INTEGER :: handle, idomain, index1_prev
INTEGER :: handle, idomain, index1_prev, sum_temp
INTEGER, ALLOCATABLE, DIMENSION(:) :: nneighbors
CALL timeset(routineN, handle)
@ -799,15 +861,19 @@ CONTAINS
WRITE (unit_nr, '(T2,A,T48,E33.3)') "eps_filter:", almo_scf_env%eps_filter
WRITE (unit_nr, '(T2,A)') "optimization of block-diagonal ALMOs:"
SELECT CASE (almo_scf_env%almo_update_algorithm)
CASE (almo_scf_diag)
! the DIIS algorith is the only choice for the diagonlaization-based algorithm
CALL print_optimizer_options(almo_scf_env%opt_block_diag_diis, unit_nr)
CASE (almo_scf_pcg)
! print out PCG options
CALL print_optimizer_options(almo_scf_env%opt_block_diag_pcg, unit_nr)
END SELECT
IF (almo_scf_env%almo_update_algorithm .EQ. almo_scf_skip) THEN
WRITE (unit_nr, '(T2,A)') "skip optimization of block-diagonal ALMOs"
ELSE
WRITE (unit_nr, '(T2,A)') "optimization of block-diagonal ALMOs:"
SELECT CASE (almo_scf_env%almo_update_algorithm)
CASE (almo_scf_diag)
! the DIIS algorith is the only choice for the diagonlaization-based algorithm
CALL print_optimizer_options(almo_scf_env%opt_block_diag_diis, unit_nr)
CASE (almo_scf_pcg)
! print out PCG options
CALL print_optimizer_options(almo_scf_env%opt_block_diag_pcg, unit_nr)
END SELECT
ENDIF
SELECT CASE (almo_scf_env%deloc_method)
CASE (almo_deloc_none)
@ -883,35 +949,118 @@ CONTAINS
! WRITE(unit_nr,'(T2,A,T48,A33)') "Parallel distribution for MOs","MOLECULAR"
!END SELECT
! print fragment's statistics
WRITE (unit_nr, '(T2,A)') REPEAT("-", 79)
WRITE (unit_nr, '(T2,A,T48,I33)') "Total fragments:", &
almo_scf_env%ndomains
IF (almo_scf_env%ndomains .LE. 200) THEN
sum_temp = SUM(almo_scf_env%nbasis_of_domain(:))
WRITE (unit_nr, '(T2,A,T53,I5,F9.2,I5,I9)') &
"Basis set size per fragment (min, av, max, total):", &
MINVAL(almo_scf_env%nbasis_of_domain(:)), &
(1.0_dp*sum_temp)/almo_scf_env%ndomains, &
MAXVAL(almo_scf_env%nbasis_of_domain(:)), &
sum_temp
!WRITE (unit_nr, '(T2,I13,F13.3,I13,I13)') &
! MINVAL(almo_scf_env%nbasis_of_domain(:)), &
! (1.0_dp*sum_temp) / almo_scf_env%ndomains, &
! MAXVAL(almo_scf_env%nbasis_of_domain(:)), &
! sum_temp
sum_temp = SUM(almo_scf_env%nocc_of_domain(:, :))
WRITE (unit_nr, '(T2,A,T53,I5,F9.2,I5,I9)') &
"Occupied MOs per fragment (min, av, max, total):", &
MINVAL(SUM(almo_scf_env%nocc_of_domain, DIM=2)), &
(1.0_dp*sum_temp)/almo_scf_env%ndomains, &
MAXVAL(SUM(almo_scf_env%nocc_of_domain, DIM=2)), &
sum_temp
!WRITE (unit_nr, '(T2,I13,F13.3,I13,I13)') &
! MINVAL( SUM(almo_scf_env%nocc_of_domain, DIM=2) ), &
! (1.0_dp*sum_temp) / almo_scf_env%ndomains, &
! MAXVAL( SUM(almo_scf_env%nocc_of_domain, DIM=2) ), &
! sum_temp
sum_temp = SUM(almo_scf_env%nvirt_of_domain(:, :))
WRITE (unit_nr, '(T2,A,T53,I5,F9.2,I5,I9)') &
"Virtual MOs per fragment (min, av, max, total):", &
MINVAL(SUM(almo_scf_env%nvirt_of_domain, DIM=2)), &
(1.0_dp*sum_temp)/almo_scf_env%ndomains, &
MAXVAL(SUM(almo_scf_env%nvirt_of_domain, DIM=2)), &
sum_temp
!WRITE (unit_nr, '(T2,I13,F13.3,I13,I13)') &
! MINVAL( SUM(almo_scf_env%nvirt_of_domain, DIM=2) ), &
! (1.0_dp*sum_temp) / almo_scf_env%ndomains, &
! MAXVAL( SUM(almo_scf_env%nvirt_of_domain, DIM=2) ), &
! sum_temp
sum_temp = SUM(almo_scf_env%charge_of_domain(:))
WRITE (unit_nr, '(T2,A,T53,I5,F9.2,I5,I9)') &
"Charges per fragment (min, av, max, total):", &
MINVAL(almo_scf_env%charge_of_domain(:)), &
(1.0_dp*sum_temp)/almo_scf_env%ndomains, &
MAXVAL(almo_scf_env%charge_of_domain(:)), &
sum_temp
!WRITE (unit_nr, '(T2,I13,F13.3,I13,I13)') &
! MINVAL(almo_scf_env%charge_of_domain(:)), &
! (1.0_dp*sum_temp) / almo_scf_env%ndomains, &
! MAXVAL(almo_scf_env%charge_of_domain(:)), &
! sum_temp
! compute the number of neighbors of each fragment
ALLOCATE (nneighbors(almo_scf_env%ndomains))
DO idomain = 1, almo_scf_env%ndomains
IF (idomain .EQ. 1) THEN
index1_prev = 1
ELSE
index1_prev = almo_scf_env%domain_map(1)%index1(idomain-1)
ENDIF
SELECT CASE (almo_scf_env%deloc_method)
CASE (almo_deloc_none)
nneighbors(idomain) = 0
CASE (almo_deloc_x, almo_deloc_scf, almo_deloc_x_then_scf)
nneighbors(idomain) = almo_scf_env%ndomains-1 ! minus self
CASE (almo_deloc_xalmo_1diag, almo_deloc_xalmo_x, almo_deloc_xalmo_scf)
nneighbors(idomain) = almo_scf_env%domain_map(1)%index1(idomain)-index1_prev-1 ! minus self
CASE DEFAULT
nneighbors(idomain) = -1
END SELECT
ENDDO ! cycle over domains
sum_temp = SUM(nneighbors(:))
WRITE (unit_nr, '(T2,A,T53,I5,F9.2,I5,I9)') &
"Deloc. neighbors of fragment (min, av, max, total):", &
MINVAL(nneighbors(:)), &
(1.0_dp*sum_temp)/almo_scf_env%ndomains, &
MAXVAL(nneighbors(:)), &
sum_temp
WRITE (unit_nr, '(T2,A)') REPEAT("-", 79)
WRITE (unit_nr, '()')
IF (almo_scf_env%ndomains .LE. 64) THEN
! print fragment info
WRITE (unit_nr, '(T2,A13,A13,A13,A13,A13,A13)') &
WRITE (unit_nr, '(T2,A10,A13,A13,A13,A13,A13)') &
"Fragment", "Basis Set", "Occupied", "Virtual", "Charge", "Deloc Neig" !,"Discarded Virt"
WRITE (unit_nr, '(T2,A)') REPEAT("-", 79)
DO idomain = 1, almo_scf_env%ndomains
IF (idomain .EQ. 1) THEN
index1_prev = 1
ELSE
index1_prev = almo_scf_env%domain_map(1)%index1(idomain-1)
ENDIF
SELECT CASE (almo_scf_env%deloc_method)
CASE (almo_deloc_none)
neig_string = "NONE"
CASE (almo_deloc_x, almo_deloc_scf, almo_deloc_x_then_scf)
neig_string = "ALL"
CASE (almo_deloc_xalmo_1diag, almo_deloc_xalmo_x, almo_deloc_xalmo_scf)
WRITE (neig_string, '(I13)') &
almo_scf_env%domain_map(1)%index1(idomain)-index1_prev-1 ! minus self
WRITE (neig_string, '(I13)') nneighbors(idomain)
CASE DEFAULT
neig_string = "N/A"
END SELECT
WRITE (unit_nr, '(T2,I13,I13,I13,I13,I13,A13)') &
WRITE (unit_nr, '(T2,I10,I13,I13,I13,I13,A13)') &
idomain, almo_scf_env%nbasis_of_domain(idomain), &
SUM(almo_scf_env%nocc_of_domain(idomain, :)), &
SUM(almo_scf_env%nvirt_of_domain(idomain, :)), &
@ -938,8 +1087,8 @@ CONTAINS
index1_prev = almo_scf_env%domain_map(1)%index1(idomain-1)
ENDIF
WRITE (unit_nr, '(T2,I13,":")') idomain
WRITE (unit_nr, '(T15,5I13)') &
WRITE (unit_nr, '(T2,I10,":")') idomain
WRITE (unit_nr, '(T12,11I6)') &
almo_scf_env%domain_map(1)%pairs &
(index1_prev:almo_scf_env%domain_map(1)%index1(idomain)-1, 1) ! includes self
@ -947,13 +1096,19 @@ CONTAINS
END SELECT
WRITE (unit_nr, '(T2,A)') REPEAT("-", 79)
ELSE ! too big to print details for each fragment
ENDIF
WRITE (unit_nr, '(T2,A)') "The system is too big to print details for each fragment."
ENDIF ! how many fragments?
WRITE (unit_nr, '(T2,A)') REPEAT("-", 79)
WRITE (unit_nr, '()')
ENDIF
DEALLOCATE (nneighbors)
ENDIF ! unit_nr > 0
CALL timestop(handle)
@ -997,9 +1152,9 @@ CONTAINS
CALL dbcsr_add_on_diag(almo_scf_env%matrix_s_blk(1), 1.0_dp)
ELSE
CALL matrix_qs_to_almo(matrix_s, almo_scf_env%matrix_s(1), &
almo_scf_env, .FALSE.)
CALL matrix_qs_to_almo(matrix_s, almo_scf_env%matrix_s_blk(1), &
almo_scf_env, .TRUE.)
almo_scf_env%mat_distr_aos, .FALSE.)
CALL dbcsr_copy(almo_scf_env%matrix_s_blk(1), &
almo_scf_env%matrix_s(1), keep_sparsity=.TRUE.)
ENDIF
CALL dbcsr_filter(almo_scf_env%matrix_s(1), almo_scf_env%eps_filter)
@ -1011,12 +1166,14 @@ CONTAINS
almo_scf_env%matrix_s_blk(1), &
threshold=almo_scf_env%eps_filter, &
order=almo_scf_env%order_lanczos, &
!order=0, &
eps_lanczos=almo_scf_env%eps_lanczos, &
max_iter_lanczos=almo_scf_env%max_iter_lanczos)
ELSE IF (almo_scf_env%almo_update_algorithm .EQ. almo_scf_dm_sign) THEN
CALL invert_Hotelling(almo_scf_env%matrix_s_blk_inv(1), &
almo_scf_env%matrix_s_blk(1), &
threshold=almo_scf_env%eps_filter)
threshold=almo_scf_env%eps_filter, &
filter_eps=almo_scf_env%eps_filter)
ENDIF
CALL timestop(handle)
@ -1051,7 +1208,8 @@ CONTAINS
unit_nr = -1
ENDIF
IF (almo_scf_env%almo_update_algorithm .EQ. almo_scf_pcg) THEN
SELECT CASE (almo_scf_env%almo_update_algorithm)
CASE (almo_scf_pcg)
! ALMO PCG optimizer as a special case of XALMO PCG
CALL almo_scf_xalmo_pcg(qs_env=qs_env, &
@ -1064,15 +1222,43 @@ CONTAINS
perturbation_only=.FALSE., &
special_case=xalmo_case_block_diag)
CALL almo_scf_t_blk_to_t_blk_orthonormal(almo_scf_env)
DO ispin = 1, almo_scf_env%nspins
CALL orthogonalize_mos(ket=almo_scf_env%matrix_t_blk(ispin), &
overlap=almo_scf_env%matrix_sigma_blk(ispin), &
metric=almo_scf_env%matrix_s_blk(1), &
retain_locality=.TRUE., &
only_normalize=.FALSE., &
nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), &
eps_filter=almo_scf_env%eps_filter, &
order_lanczos=almo_scf_env%order_lanczos, &
eps_lanczos=almo_scf_env%eps_lanczos, &
max_iter_lanczos=almo_scf_env%max_iter_lanczos)
ENDDO
!CALL almo_scf_t_blk_to_t_blk_orthonormal(almo_scf_env)
ELSE
CASE (almo_scf_diag)
! mixing/DIIS optimizer
CALL almo_scf_block_diagonal(qs_env, almo_scf_env, &
almo_scf_env%opt_block_diag_diis)
ENDIF
CASE (almo_scf_skip)
DO ispin = 1, almo_scf_env%nspins
CALL orthogonalize_mos(ket=almo_scf_env%matrix_t_blk(ispin), &
overlap=almo_scf_env%matrix_sigma_blk(ispin), &
metric=almo_scf_env%matrix_s_blk(1), &
retain_locality=.TRUE., &
only_normalize=.FALSE., &
nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), &
eps_filter=almo_scf_env%eps_filter, &
order_lanczos=almo_scf_env%order_lanczos, &
eps_lanczos=almo_scf_env%eps_lanczos, &
max_iter_lanczos=almo_scf_env%max_iter_lanczos)
ENDDO
!CALL almo_scf_t_blk_to_t_blk_orthonormal(almo_scf_env)
END SELECT
! we might need a copy of the converged KS and sigma_inv
DO ispin = 1, almo_scf_env%nspins
@ -1220,7 +1406,8 @@ CONTAINS
quench_t=no_quench, &
matrix_t_in=almo_scf_env%matrix_t_blk, &
matrix_t_out=almo_scf_env%matrix_t, &
assume_t0_q0x=.TRUE., &
!assume_t0_q0x=.TRUE., &
assume_t0_q0x=(almo_scf_env%xalmo_trial_wf .EQ. xalmo_trial_r0_out), &
perturbation_only=.TRUE., &
special_case=xalmo_case_fully_deloc)
@ -1252,7 +1439,8 @@ CONTAINS
quench_t=almo_scf_env%quench_t, &
matrix_t_in=almo_scf_env%matrix_t_blk, &
matrix_t_out=almo_scf_env%matrix_t, &
assume_t0_q0x=.TRUE., &
!assume_t0_q0x=.TRUE., &
assume_t0_q0x=(almo_scf_env%xalmo_trial_wf .EQ. xalmo_trial_r0_out), &
perturbation_only=.TRUE., &
special_case=xalmo_case_normal)
@ -1284,7 +1472,8 @@ CONTAINS
quench_t=almo_scf_env%quench_t, &
matrix_t_in=almo_scf_env%matrix_t_blk, &
matrix_t_out=almo_scf_env%matrix_t, &
assume_t0_q0x=.TRUE., &
!assume_t0_q0x=.TRUE., &
assume_t0_q0x=(almo_scf_env%xalmo_trial_wf .EQ. xalmo_trial_r0_out), &
perturbation_only=.FALSE., &
special_case=xalmo_case_normal)
@ -1358,39 +1547,44 @@ CONTAINS
INTEGER :: handle, ispin
TYPE(cp_fm_type), POINTER :: mo_coeff
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_w
TYPE(dbcsr_type), ALLOCATABLE, DIMENSION(:) :: matrix_t_orthogonal
TYPE(dbcsr_type), ALLOCATABLE, DIMENSION(:) :: matrix_t_processed
TYPE(mo_set_p_type), DIMENSION(:), POINTER :: mos
TYPE(qs_scf_env_type), POINTER :: scf_env
CALL timeset(routineN, handle)
! store the matrices for a next scf run
CALL almo_scf_store_result(almo_scf_env)
! store matrices to speed up the next scf run
CALL almo_scf_store_extrapolation_data(almo_scf_env)
! orthogonalize orbitals before returning them to QS
ALLOCATE (matrix_t_orthogonal(almo_scf_env%nspins))
ALLOCATE (matrix_t_processed(almo_scf_env%nspins))
DO ispin = 1, almo_scf_env%nspins
CALL dbcsr_create(matrix_t_orthogonal(ispin), &
CALL dbcsr_create(matrix_t_processed(ispin), &
template=almo_scf_env%matrix_t(ispin), &
matrix_type=dbcsr_type_no_symmetry)
CALL dbcsr_copy(matrix_t_orthogonal(ispin), &
CALL dbcsr_copy(matrix_t_processed(ispin), &
almo_scf_env%matrix_t(ispin))
CALL orthogonalize_mos(ket=matrix_t_orthogonal(ispin), &
overlap=almo_scf_env%matrix_sigma(ispin), &
metric=almo_scf_env%matrix_s(1), &
retain_locality=.FALSE., &
only_normalize=.FALSE., &
eps_filter=almo_scf_env%eps_filter, &
order_lanczos=almo_scf_env%order_lanczos, &
eps_lanczos=almo_scf_env%eps_lanczos, &
max_iter_lanczos=almo_scf_env%max_iter_lanczos)
IF (almo_scf_env%return_orthogonalized_mos) THEN
CALL orthogonalize_mos(ket=matrix_t_processed(ispin), &
overlap=almo_scf_env%matrix_sigma(ispin), &
metric=almo_scf_env%matrix_s(1), &
retain_locality=.FALSE., &
only_normalize=.FALSE., &
nocc_of_domain=almo_scf_env%nocc_of_domain(:, ispin), &
eps_filter=almo_scf_env%eps_filter, &
order_lanczos=almo_scf_env%order_lanczos, &
eps_lanczos=almo_scf_env%eps_lanczos, &
max_iter_lanczos=almo_scf_env%max_iter_lanczos)
ENDIF
ENDDO
! return orthogonalized orbitals to QS
! return orbitals to QS
NULLIFY (mos, mo_coeff, scf_env)
CALL get_qs_env(qs_env, mos=mos, scf_env=scf_env)
@ -1399,12 +1593,10 @@ CONTAINS
! Currently only fm version of mo_set is usable.
! First transform the matrix_t to fm version
CALL get_mo_set(mos(ispin)%mo_set, mo_coeff=mo_coeff)
CALL copy_dbcsr_to_fm(matrix_t_orthogonal(ispin), mo_coeff)
CALL dbcsr_release(matrix_t_orthogonal(ispin))
CALL copy_dbcsr_to_fm(matrix_t_processed(ispin), mo_coeff)
CALL dbcsr_release(matrix_t_processed(ispin))
ENDDO
DEALLOCATE (matrix_t_orthogonal)
DEALLOCATE (matrix_t_processed)
! calculate post scf properties
! CALL almo_post_scf_compute_properties(qs_env, almo_scf_env)
@ -2093,7 +2285,6 @@ CONTAINS
!> \author Yifei Shi
! **************************************************************************************************
SUBROUTINE almo_post_scf_compute_properties(qs_env)
! SUBROUTINE almo_post_scf_compute_properties(qs_env, almo_scf_env)
TYPE(qs_environment_type), POINTER :: qs_env
CHARACTER(len=*), PARAMETER :: routineN = 'almo_post_scf_compute_properties', &
@ -2102,8 +2293,6 @@ CONTAINS
TYPE(qs_scf_env_type), POINTER :: scf_env
TYPE(section_vals_type), POINTER :: dft_section, input
!TYPE(almo_scf_env_type) :: almo_scf_env
CALL get_qs_env(qs_env, scf_env=scf_env, input=input)
dft_section => section_vals_get_subs_vals(input, "DFT")
CALL qs_scf_compute_properties(qs_env, dft_section)

91
src/almo_scf_env_methods.F
View file

@ -15,10 +15,12 @@ MODULE almo_scf_env_methods
almo_scf_env_type
USE cp_control_types, ONLY: dft_control_type
USE input_constants, ONLY: &
almo_constraint_distance, almo_deloc_xalmo_1diag, almo_domain_layout_atomic, &
almo_domain_layout_molecular, almo_frz_crystal, almo_mat_distr_molecular, almo_scf_diag, &
almo_scf_pcg, cg_hager_zhang, do_bondparm_vdw, molecular_guess, tensor_orthogonal, &
virt_full, virt_minimal, virt_number
almo_constraint_distance, almo_deloc_none, almo_deloc_scf, almo_deloc_x, &
almo_deloc_x_then_scf, almo_deloc_xalmo_1diag, almo_domain_layout_atomic, &
almo_domain_layout_molecular, almo_frz_crystal, almo_mat_distr_molecular, &
almo_occ_vol_penalty_none, almo_scf_diag, almo_scf_pcg, almo_scf_skip, cg_hager_zhang, &
do_bondparm_vdw, molecular_guess, tensor_orthogonal, virt_full, virt_minimal, virt_number, &
xalmo_trial_r0_out
USE input_section_types, ONLY: section_vals_get_subs_vals,&
section_vals_type,&
section_vals_val_get
@ -105,7 +107,8 @@ CONTAINS
INTEGER :: handle
TYPE(section_vals_type), POINTER :: almo_analysis_section, almo_opt_diis_section, &
almo_opt_pcg_section, almo_scf_section, xalmo_opt_pcg_section
almo_opt_pcg_section, almo_scf_section, matrix_iterate_section, penalty_section, &
xalmo_opt_newton_pcg_section, xalmo_opt_pcg_section
CALL timeset(routineN, handle)
@ -117,12 +120,11 @@ CONTAINS
xalmo_opt_pcg_section => section_vals_get_subs_vals(almo_scf_section, &
"XALMO_OPTIMIZER_PCG")
almo_analysis_section => section_vals_get_subs_vals(almo_scf_section, "ANALYSIS")
! RZK-warning the values of these keywords are hardcoded
! but ideally they should also be read from the input file
almo_scf_env%order_lanczos = 3
almo_scf_env%eps_lanczos = 1.0E-4_dp
almo_scf_env%max_iter_lanczos = 40
xalmo_opt_newton_pcg_section => section_vals_get_subs_vals(xalmo_opt_pcg_section, &
"XALMO_NEWTON_PCG_SOLVER")
matrix_iterate_section => section_vals_get_subs_vals(almo_scf_section, &
"MATRIX_ITERATE")
penalty_section => section_vals_get_subs_vals(almo_scf_section, "PENALTY")
! read user input
! common ALMO options
@ -132,6 +134,10 @@ CONTAINS
i_val=almo_scf_env%almo_scf_guess)
CALL section_vals_val_get(almo_scf_section, "ALMO_ALGORITHM", &
i_val=almo_scf_env%almo_update_algorithm)
CALL section_vals_val_get(almo_scf_section, "XALMO_TRIAL_WF", &
i_val=almo_scf_env%xalmo_trial_wf)
CALL section_vals_val_get(almo_scf_section, "MO_OVERLAP_INV_ALG", &
i_val=almo_scf_env%sigma_inv_algorithm)
CALL section_vals_val_get(almo_scf_section, "DELOCALIZE_METHOD", &
i_val=almo_scf_env%deloc_method)
CALL section_vals_val_get(almo_scf_section, "XALMO_R_CUTOFF_FACTOR", &
@ -142,12 +148,27 @@ CONTAINS
CALL section_vals_val_get(almo_scf_section, "XALMO_EXTRAPOLATION_ORDER", &
i_val=almo_scf_env%xalmo_extrapolation_order)
almo_scf_env%xalmo_extrapolation_order = MAX(0, almo_scf_env%xalmo_extrapolation_order)
CALL section_vals_val_get(almo_scf_section, "RETURN_ORTHOGONALIZED_MOS", &
l_val=almo_scf_env%return_orthogonalized_mos)
CALL section_vals_val_get(matrix_iterate_section, "EPS_LANCZOS", &
r_val=almo_scf_env%eps_lanczos)
CALL section_vals_val_get(matrix_iterate_section, "ORDER_LANCZOS", &
i_val=almo_scf_env%order_lanczos)
CALL section_vals_val_get(matrix_iterate_section, "MAX_ITER_LANCZOS", &
i_val=almo_scf_env%max_iter_lanczos)
CALL section_vals_val_get(matrix_iterate_section, "EPS_TARGET_FACTOR", &
r_val=almo_scf_env%matrix_iter_eps_error_factor)
! optimizers
CALL section_vals_val_get(almo_opt_diis_section, "EPS_ERROR", &
r_val=almo_scf_env%opt_block_diag_diis%eps_error)
CALL section_vals_val_get(almo_opt_diis_section, "MAX_ITER", &
i_val=almo_scf_env%opt_block_diag_diis%max_iter)
CALL section_vals_val_get(almo_opt_diis_section, "EPS_ERROR_EARLY", &
r_val=almo_scf_env%opt_block_diag_diis%eps_error_early)
CALL section_vals_val_get(almo_opt_diis_section, "MAX_ITER_EARLY", &
i_val=almo_scf_env%opt_block_diag_diis%max_iter_early)
CALL section_vals_val_get(almo_opt_diis_section, "N_DIIS", &
i_val=almo_scf_env%opt_block_diag_diis%ndiis)
@ -155,6 +176,10 @@ CONTAINS
r_val=almo_scf_env%opt_block_diag_pcg%eps_error)
CALL section_vals_val_get(almo_opt_pcg_section, "MAX_ITER", &
i_val=almo_scf_env%opt_block_diag_pcg%max_iter)
CALL section_vals_val_get(almo_opt_pcg_section, "EPS_ERROR_EARLY", &
r_val=almo_scf_env%opt_block_diag_pcg%eps_error_early)
CALL section_vals_val_get(almo_opt_pcg_section, "MAX_ITER_EARLY", &
i_val=almo_scf_env%opt_block_diag_pcg%max_iter_early)
CALL section_vals_val_get(almo_opt_pcg_section, "MAX_ITER_OUTER_LOOP", &
i_val=almo_scf_env%opt_block_diag_pcg%max_iter_outer_loop)
CALL section_vals_val_get(almo_opt_pcg_section, "LIN_SEARCH_EPS_ERROR", &
@ -170,17 +195,39 @@ CONTAINS
r_val=almo_scf_env%opt_xalmo_pcg%eps_error)
CALL section_vals_val_get(xalmo_opt_pcg_section, "MAX_ITER", &
i_val=almo_scf_env%opt_xalmo_pcg%max_iter)
CALL section_vals_val_get(xalmo_opt_pcg_section, "EPS_ERROR_EARLY", &
r_val=almo_scf_env%opt_xalmo_pcg%eps_error_early)
CALL section_vals_val_get(xalmo_opt_pcg_section, "MAX_ITER_EARLY", &
i_val=almo_scf_env%opt_xalmo_pcg%max_iter_early)
CALL section_vals_val_get(xalmo_opt_pcg_section, "MAX_ITER_OUTER_LOOP", &
i_val=almo_scf_env%opt_xalmo_pcg%max_iter_outer_loop)
CALL section_vals_val_get(xalmo_opt_pcg_section, "LIN_SEARCH_EPS_ERROR", &
r_val=almo_scf_env%opt_xalmo_pcg%lin_search_eps_error)
CALL section_vals_val_get(xalmo_opt_pcg_section, "LIN_SEARCH_STEP_SIZE_GUESS", &
r_val=almo_scf_env%opt_xalmo_pcg%lin_search_step_size_guess)
CALL section_vals_val_get(xalmo_opt_pcg_section, "PRECOND_FILTER_THRESHOLD", &
r_val=almo_scf_env%opt_xalmo_pcg%neglect_threshold)
CALL section_vals_val_get(xalmo_opt_pcg_section, "CONJUGATOR", &
i_val=almo_scf_env%opt_xalmo_pcg%conjugator)
CALL section_vals_val_get(xalmo_opt_pcg_section, "PRECONDITIONER", &
i_val=almo_scf_env%opt_xalmo_pcg%preconditioner)
CALL section_vals_val_get(xalmo_opt_newton_pcg_section, "EPS_ERROR", &
r_val=almo_scf_env%opt_xalmo_newton_pcg_solver%eps_error)
CALL section_vals_val_get(xalmo_opt_newton_pcg_section, "MAX_ITER", &
i_val=almo_scf_env%opt_xalmo_newton_pcg_solver%max_iter)
CALL section_vals_val_get(xalmo_opt_newton_pcg_section, "MAX_ITER_OUTER_LOOP", &
i_val=almo_scf_env%opt_xalmo_newton_pcg_solver%max_iter_outer_loop)
CALL section_vals_val_get(xalmo_opt_newton_pcg_section, "PRECONDITIONER", &
i_val=almo_scf_env%opt_xalmo_newton_pcg_solver%preconditioner)
CALL section_vals_val_get(penalty_section, &
"OCCUPIED_VOLUME_PENALTY_COEFF", &
r_val=almo_scf_env%penalty%occ_vol_coeff)
CALL section_vals_val_get(penalty_section, &
"OCCUPIED_VOLUME_PENALTY_METHOD", &
i_val=almo_scf_env%penalty%occ_vol_method)
CALL section_vals_val_get(almo_analysis_section, "_SECTION_PARAMETERS_", &
l_val=almo_scf_env%almo_analysis%do_analysis)
CALL section_vals_val_get(almo_analysis_section, "FROZEN_MO_ENERGY_TERM", &
@ -296,8 +343,8 @@ CONTAINS
almo_scf_env%constraint_type = almo_constraint_distance
almo_scf_env%mu = -0.1_dp
almo_scf_env%fixed_mu = .FALSE.
almo_scf_env%eps_prev_guess = almo_scf_env%eps_filter/100.0_dp
almo_scf_env%mixing_fraction = 0.45_dp
almo_scf_env%eps_prev_guess = almo_scf_env%eps_filter/1000.0_dp
almo_scf_env%deloc_cayley_tensor_type = tensor_orthogonal
almo_scf_env%deloc_cayley_conjugator = cg_hager_zhang
@ -362,6 +409,17 @@ CONTAINS
ENDIF
! check for conflicts between options
IF (almo_scf_env%xalmo_trial_wf .EQ. xalmo_trial_r0_out .AND. &
almo_scf_env%almo_update_algorithm .EQ. almo_scf_skip .AND. &
almo_scf_env%almo_scf_guess .NE. molecular_guess) THEN
CPABORT("R0 projector requires optimized ALMOs")
ENDIF
IF (almo_scf_env%deloc_method .EQ. almo_deloc_none .AND. &
almo_scf_env%almo_update_algorithm .EQ. almo_scf_skip) THEN
CPABORT("No optimization requested")
ENDIF
IF (almo_scf_env%deloc_truncate_virt .EQ. virt_number .AND. &
almo_scf_env%deloc_virt_per_domain .LE. 0) THEN
CPABORT("specify a positive number of virtual orbitals")
@ -409,6 +467,15 @@ CONTAINS
ENDIF ! end analysis settings
! check penalty settings
IF (almo_scf_env%penalty%occ_vol_method .NE. almo_occ_vol_penalty_none) THEN
IF (almo_scf_env%deloc_method .NE. almo_deloc_x .AND. &
almo_scf_env%deloc_method .NE. almo_deloc_scf .AND. &
almo_scf_env%deloc_method .NE. almo_deloc_x_then_scf) THEN
CPABORT("Occupied volume penalty seems to work only with completely delocalized orbitals")
ENDIF
ENDIF ! end penalty settings
CALL timestop(handle)
END SUBROUTINE almo_scf_init_read_write_input

1721
src/almo_scf_methods.F
View file

File diff suppressed because it is too large Load diff

4555
src/almo_scf_optimizer.F
View file

File diff suppressed because it is too large Load diff

311
src/almo_scf_qs.F
View file

@ -72,7 +72,6 @@ MODULE almo_scf_qs
USE qs_rho_methods, ONLY: qs_rho_update_rho
USE qs_rho_types, ONLY: qs_rho_get,&
qs_rho_type
USE qs_scf_initialization, ONLY: qs_scf_env_initialize
USE qs_scf_types, ONLY: qs_scf_env_type,&
scf_env_create
#include "./base/base_uses.f90"
@ -83,10 +82,15 @@ MODULE almo_scf_qs
CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'almo_scf_qs'
PUBLIC :: almo_scf_init_qs, matrix_almo_create, matrix_qs_to_almo, &
almo_scf_dm_to_ks, calculate_w_matrix_almo, &
almo_scf_construct_quencher, almo_scf_update_ks_energy, &
matrix_almo_to_qs
PUBLIC :: matrix_almo_create, &
almo_scf_construct_quencher, &
calculate_w_matrix_almo, &
init_almo_ks_matrix_via_qs, &
almo_scf_update_ks_energy, &
construct_qs_mos, &
matrix_qs_to_almo, &
almo_dm_to_almo_ks, &
almo_dm_to_qs_env
CONTAINS
@ -125,7 +129,7 @@ CONTAINS
INTEGER, DIMENSION(:), POINTER :: blk_distr, blk_sizes, block_sizes_new, cluster_distr, &
cluster_distr_new, col_cluster_new, col_distr_new, col_sizes_new, distr_new_array, &
row_cluster_new, row_distr_new, row_sizes_new
LOGICAL :: active, tr
LOGICAL :: active, one_dim_is_mo, tr
REAL(KIND=dp), DIMENSION(:, :), POINTER :: p_new_block
TYPE(dbcsr_distribution_type) :: dist_new, dist_qs
@ -233,6 +237,9 @@ CONTAINS
ALLOCATE (block_sizes_new(nlength))
IF (size_keys(dimen) == almo_mat_dim_occ) THEN
block_sizes_new(:) = almo_scf_env%nocc_of_domain(:, spin_key)
! Handle zero-electron fragments by adding one-orbital that
! must remain zero at all times
WHERE (block_sizes_new == 0) block_sizes_new = 1
ELSE IF (size_keys(dimen) == almo_mat_dim_virt_disc) THEN
block_sizes_new(:) = almo_scf_env%nvirt_disc_of_domain(:, spin_key)
ELSE IF (size_keys(dimen) == almo_mat_dim_virt_full) THEN
@ -363,6 +370,7 @@ CONTAINS
!QQQ ENDIF ! mynode
!QQQ ENDDO
!QQQENDDO
!QQQtake care of zero-electron fragments
! endQQQ - end of the quadratic part
! start linear-scaling replacement:
! works only for molecular blocks AND molecular distributions
@ -377,6 +385,30 @@ CONTAINS
active = .TRUE.
one_dim_is_mo = .FALSE.
DO dimen = 1, 2 ! 1 - row, 2 - column dimension
IF (size_keys(dimen) == almo_mat_dim_occ) one_dim_is_mo = .TRUE.
ENDDO
IF (one_dim_is_mo .AND. almo_scf_env%nocc_of_domain(row, spin_key) == 0) active = .FALSE.
one_dim_is_mo = .FALSE.
DO dimen = 1, 2
IF (size_keys(dimen) == almo_mat_dim_virt) one_dim_is_mo = .TRUE.
ENDDO
IF (one_dim_is_mo .AND. almo_scf_env%nvirt_of_domain(row, spin_key) == 0) active = .FALSE.
one_dim_is_mo = .FALSE.
DO dimen = 1, 2
IF (size_keys(dimen) == almo_mat_dim_virt_disc) one_dim_is_mo = .TRUE.
ENDDO
IF (one_dim_is_mo .AND. almo_scf_env%nvirt_disc_of_domain(row, spin_key) == 0) active = .FALSE.
one_dim_is_mo = .FALSE.
DO dimen = 1, 2
IF (size_keys(dimen) == almo_mat_dim_virt_full) one_dim_is_mo = .TRUE.
ENDDO
IF (one_dim_is_mo .AND. almo_scf_env%nvirt_full_of_domain(row, spin_key) == 0) active = .FALSE.
IF (active) THEN
NULLIFY (p_new_block)
CALL dbcsr_reserve_block2d(matrix_new, iblock_row, iblock_col, p_new_block)
@ -400,16 +432,16 @@ CONTAINS
!> \brief convert between two types of matrices: QS style to ALMO style
!> \param matrix_qs ...
!> \param matrix_almo ...
!> \param almo_scf_env ...
!> \param mat_distr_aos ...
!> \param keep_sparsity ...
!> \par History
!> 2011.06 created [Rustam Z Khaliullin]
!> \author Rustam Z Khaliullin
! **************************************************************************************************
SUBROUTINE matrix_qs_to_almo(matrix_qs, matrix_almo, almo_scf_env, &
keep_sparsity)
SUBROUTINE matrix_qs_to_almo(matrix_qs, matrix_almo, mat_distr_aos, keep_sparsity)
TYPE(dbcsr_type) :: matrix_qs, matrix_almo
TYPE(almo_scf_env_type) :: almo_scf_env
INTEGER :: mat_distr_aos
LOGICAL, INTENT(IN) :: keep_sparsity
CHARACTER(len=*), PARAMETER :: routineN = 'matrix_qs_to_almo', &
@ -421,7 +453,7 @@ CONTAINS
CALL timeset(routineN, handle)
!RZK-warning if it's not a N(AO)xN(AO) matrix then stop
SELECT CASE (almo_scf_env%mat_distr_aos)
SELECT CASE (mat_distr_aos)
CASE (almo_mat_distr_atomic)
! automatic data_type conversion
CALL dbcsr_copy(matrix_almo, matrix_qs, &
@ -456,14 +488,14 @@ CONTAINS
!> \brief convert between two types of matrices: ALMO style to QS style
!> \param matrix_almo ...
!> \param matrix_qs ...
!> \param almo_scf_env ...
!> \param mat_distr_aos ...
!> \par History
!> 2011.06 created [Rustam Z Khaliullin]
!> \author Rustam Z Khaliullin
! **************************************************************************************************
SUBROUTINE matrix_almo_to_qs(matrix_almo, matrix_qs, almo_scf_env)
SUBROUTINE matrix_almo_to_qs(matrix_almo, matrix_qs, mat_distr_aos)
TYPE(dbcsr_type) :: matrix_almo, matrix_qs
TYPE(almo_scf_env_type), INTENT(IN) :: almo_scf_env
INTEGER, INTENT(IN) :: mat_distr_aos
CHARACTER(len=*), PARAMETER :: routineN = 'matrix_almo_to_qs', &
routineP = moduleN//':'//routineN
@ -473,66 +505,46 @@ CONTAINS
CALL timeset(routineN, handle)
! RZK-warning if it's not a N(AO)xN(AO) matrix then stop
! IF (ls_mstruct%single_precision) THEN
! CALL dbcsr_init (matrix_tmp)
! CALL dbcsr_create (matrix_tmp, template=matrix_ls,&
! data_type=dbcsr_type_real_8)
! CALL dbcsr_copy (matrix_tmp, matrix_ls)
! ENDIF
SELECT CASE (almo_scf_env%mat_distr_aos)
SELECT CASE (mat_distr_aos)
CASE (almo_mat_distr_atomic)
! IF (ls_mstruct%single_precision) THEN
! CALL dbcsr_copy_into_existing (matrix_qs, matrix_tmp)
! ELSE
CALL dbcsr_copy_into_existing(matrix_qs, matrix_almo)
! ENDIF
CASE (almo_mat_distr_molecular)
CALL dbcsr_set(matrix_qs, 0.0_dp)
! IF (ls_mstruct%single_precision) THEN
! CALL dbcsr_complete_redistribute(matrix_tmp, matrix_qs, keep_sparsity=.TRUE.)
! ELSE
CALL dbcsr_complete_redistribute(matrix_almo, matrix_qs, keep_sparsity=.TRUE.)
! ENDIF
CASE DEFAULT
CPABORT("")
END SELECT
! IF (ls_mstruct%single_precision) THEN
! CALL dbcsr_release(matrix_tmp)
! ENDIF
CALL timestop(handle)
END SUBROUTINE matrix_almo_to_qs
! **************************************************************************************************
!> \brief Initialization of QS and ALMOs
!> Some parts can be factored-out since they are common
!> for the other SCF methods
!> \brief Initialization of the QS and ALMO KS matrix
!> \param qs_env ...
!> \param almo_scf_env ...
!> \param matrix_ks ...
!> \param mat_distr_aos ...
!> \param eps_filter ...
!> \par History
!> 2011.05 created [Rustam Z Khaliullin]
!> \author Rustam Z Khaliullin
! **************************************************************************************************
SUBROUTINE almo_scf_init_qs(qs_env, almo_scf_env)
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(almo_scf_env_type) :: almo_scf_env
SUBROUTINE init_almo_ks_matrix_via_qs(qs_env, matrix_ks, mat_distr_aos, eps_filter)
CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_init_qs', &
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(dbcsr_type), DIMENSION(:) :: matrix_ks
INTEGER :: mat_distr_aos
REAL(KIND=dp) :: eps_filter
CHARACTER(len=*), PARAMETER :: routineN = 'init_almo_ks_matrix_via_qs', &
routineP = moduleN//':'//routineN
INTEGER :: handle, ispin, ncol_fm, nrow_fm, nspin
TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
TYPE(cp_fm_type), POINTER :: mo_fm_copy
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s
INTEGER :: handle, ispin, nspin
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_qs_ks, matrix_qs_s
TYPE(dft_control_type), POINTER :: dft_control
TYPE(mo_set_p_type), DIMENSION(:), POINTER :: mos
TYPE(neighbor_list_set_p_type), DIMENSION(:), &
POINTER :: sab_orb
TYPE(qs_ks_env_type), POINTER :: ks_env
TYPE(qs_scf_env_type), POINTER :: scf_env
CALL timeset(routineN, handle)
@ -541,33 +553,69 @@ CONTAINS
! get basic quantities from the qs_env
CALL get_qs_env(qs_env, &
dft_control=dft_control, &
matrix_s=matrix_s, &
matrix_ks=matrix_ks, &
matrix_s=matrix_qs_s, &
matrix_ks=matrix_qs_ks, &
ks_env=ks_env, &
sab_orb=sab_orb)
nspin = dft_control%nspins
! create matrix_ks if necessary
IF (.NOT. ASSOCIATED(matrix_ks)) THEN
CALL dbcsr_allocate_matrix_set(matrix_ks, nspin)
! create matrix_ks in the QS env if necessary
IF (.NOT. ASSOCIATED(matrix_qs_ks)) THEN
CALL dbcsr_allocate_matrix_set(matrix_qs_ks, nspin)
DO ispin = 1, nspin
ALLOCATE (matrix_ks(ispin)%matrix)
CALL dbcsr_create(matrix_ks(ispin)%matrix, &
template=matrix_s(1)%matrix)
CALL cp_dbcsr_alloc_block_from_nbl(matrix_ks(ispin)%matrix, sab_orb)
CALL dbcsr_set(matrix_ks(ispin)%matrix, 0.0_dp)
ALLOCATE (matrix_qs_ks(ispin)%matrix)
CALL dbcsr_create(matrix_qs_ks(ispin)%matrix, &
template=matrix_qs_s(1)%matrix)
CALL cp_dbcsr_alloc_block_from_nbl(matrix_qs_ks(ispin)%matrix, sab_orb)
CALL dbcsr_set(matrix_qs_ks(ispin)%matrix, 0.0_dp)
ENDDO
CALL set_ks_env(ks_env, matrix_ks=matrix_ks)
CALL set_ks_env(ks_env, matrix_ks=matrix_qs_ks)
ENDIF
! copy to ALMO
DO ispin = 1, nspin
CALL matrix_qs_to_almo(matrix_qs_ks(ispin)%matrix, &
matrix_ks(ispin), mat_distr_aos, .FALSE.)
CALL dbcsr_filter(matrix_ks(ispin), eps_filter)
ENDDO
CALL timestop(handle)
END SUBROUTINE init_almo_ks_matrix_via_qs
! **************************************************************************************************
!> \brief Create MOs in the QS env to be able to return ALMOs to QS
!> \param qs_env ...
!> \param almo_scf_env ...
!> \par History
!> 2016.12 created [Yifei Shi]
!> \author Yifei Shi
! **************************************************************************************************
SUBROUTINE construct_qs_mos(qs_env, almo_scf_env)
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env
CHARACTER(len=*), PARAMETER :: routineN = 'construct_qs_mos', &
routineP = moduleN//':'//routineN
INTEGER :: handle, ispin, ncol_fm, nrow_fm
TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
TYPE(cp_fm_type), POINTER :: mo_fm_copy
TYPE(dft_control_type), POINTER :: dft_control
TYPE(mo_set_p_type), DIMENSION(:), POINTER :: mos
TYPE(qs_scf_env_type), POINTER :: scf_env
CALL timeset(routineN, handle)
! create and init scf_env (this is necessary to return MOs to qs)
NULLIFY (mos, mo_fm_copy, fm_struct_tmp, scf_env)
CALL scf_env_create(scf_env)
CALL qs_scf_env_initialize(qs_env, scf_env)
!CALL qs_scf_env_initialize(qs_env, scf_env)
CALL set_qs_env(qs_env, scf_env=scf_env)
CALL get_qs_env(qs_env, mos=mos)
CALL get_qs_env(qs_env, dft_control=dft_control, mos=mos)
CALL dbcsr_get_info(almo_scf_env%matrix_t(1), nfullrows_total=nrow_fm, nfullcols_total=ncol_fm)
@ -600,55 +648,131 @@ CONTAINS
CALL timestop(handle)
END SUBROUTINE almo_scf_init_qs
END SUBROUTINE construct_qs_mos
! **************************************************************************************************
!> \brief use the density matrix in almo_scf_env
!> to compute the new energy and KS matrix
!> \brief return density matrix to the qs_env
!> \param qs_env ...
!> \param almo_scf_env ...
!> \param energy_new ...
!> \param matrix_p ...
!> \param mat_distr_aos ...
!> \par History
!> 2011.05 created [Rustam Z Khaliullin]
!> \author Rustam Z Khaliullin
! **************************************************************************************************
SUBROUTINE almo_scf_dm_to_ks(qs_env, almo_scf_env, energy_new)
SUBROUTINE almo_dm_to_qs_env(qs_env, matrix_p, mat_distr_aos)
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(almo_scf_env_type) :: almo_scf_env
REAL(KIND=dp) :: energy_new
TYPE(dbcsr_type), DIMENSION(:) :: matrix_p
INTEGER, INTENT(IN) :: mat_distr_aos
CHARACTER(len=*), PARAMETER :: routineN = 'almo_scf_dm_to_ks', &
CHARACTER(len=*), PARAMETER :: routineN = 'almo_dm_to_qs_env', &
routineP = moduleN//':'//routineN
INTEGER :: handle, ispin, nspin
INTEGER :: handle, ispin, nspins
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao
TYPE(qs_energy_type), POINTER :: energy
TYPE(qs_rho_type), POINTER :: rho
NULLIFY (rho, rho_ao)
CALL timeset(routineN, handle)
nspin = almo_scf_env%nspins
CALL get_qs_env(qs_env, energy=energy, rho=rho)
NULLIFY (rho, rho_ao)
nspins = SIZE(matrix_p)
CALL get_qs_env(qs_env, rho=rho)
CALL qs_rho_get(rho, rho_ao=rho_ao)
! set the new density matrix
DO ispin = 1, nspin
CALL matrix_almo_to_qs(almo_scf_env%matrix_p(ispin), &
DO ispin = 1, nspins
CALL matrix_almo_to_qs(matrix_p(ispin), &
rho_ao(ispin)%matrix, &
almo_scf_env)
mat_distr_aos)
END DO
! compute the corresponding KS matrix and new energy
CALL qs_rho_update_rho(rho, qs_env=qs_env)
CALL qs_ks_did_change(qs_env%ks_env, rho_changed=.TRUE.)
CALL qs_ks_update_qs_env(qs_env, calculate_forces=.FALSE., just_energy=.FALSE., &
print_active=.TRUE.)
energy_new = energy%total
CALL timestop(handle)
END SUBROUTINE almo_scf_dm_to_ks
END SUBROUTINE almo_dm_to_qs_env
! **************************************************************************************************
!> \brief uses the ALMO density matrix
!> to compute KS matrix (inside QS environment) and the new energy
!> \param qs_env ...
!> \param matrix_p ...
!> \param energy_total ...
!> \param mat_distr_aos ...
!> \par History
!> 2011.05 created [Rustam Z Khaliullin]
!> \author Rustam Z Khaliullin
! **************************************************************************************************
SUBROUTINE almo_dm_to_qs_ks(qs_env, matrix_p, energy_total, mat_distr_aos)
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(dbcsr_type), DIMENSION(:) :: matrix_p
REAL(KIND=dp) :: energy_total
INTEGER, INTENT(IN) :: mat_distr_aos
CHARACTER(len=*), PARAMETER :: routineN = 'almo_dm_to_qs_ks', &
routineP = moduleN//':'//routineN
INTEGER :: handle
TYPE(qs_energy_type), POINTER :: energy
CALL timeset(routineN, handle)
NULLIFY (energy)
CALL get_qs_env(qs_env, energy=energy)
CALL almo_dm_to_qs_env(qs_env, matrix_p, mat_distr_aos)
CALL qs_ks_update_qs_env(qs_env, calculate_forces=.FALSE., just_energy=.FALSE., &
print_active=.TRUE.)
energy_total = energy%total
CALL timestop(handle)
END SUBROUTINE almo_dm_to_qs_ks
! **************************************************************************************************
!> \brief uses the ALMO density matrix
!> to compute ALMO KS matrix and the new energy
!> \param qs_env ...
!> \param matrix_p ...
!> \param matrix_ks ...
!> \param energy_total ...
!> \param eps_filter ...
!> \param mat_distr_aos ...
!> \par History
!> 2011.05 created [Rustam Z Khaliullin]
!> \author Rustam Z Khaliullin
! **************************************************************************************************
SUBROUTINE almo_dm_to_almo_ks(qs_env, matrix_p, matrix_ks, energy_total, eps_filter, &
mat_distr_aos)
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(dbcsr_type), DIMENSION(:) :: matrix_p, matrix_ks
REAL(KIND=dp) :: energy_total, eps_filter
INTEGER, INTENT(IN) :: mat_distr_aos
CHARACTER(len=*), PARAMETER :: routineN = 'almo_dm_to_almo_ks', &
routineP = moduleN//':'//routineN
INTEGER :: handle, ispin, nspins
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_qs_ks
CALL timeset(routineN, handle)
! update KS matrix in the QS env
CALL almo_dm_to_qs_ks(qs_env, matrix_p, energy_total, mat_distr_aos)
nspins = SIZE(matrix_ks)
! get KS matrix from the QS env and convert to the ALMO format
CALL get_qs_env(qs_env, matrix_ks=matrix_qs_ks)
DO ispin = 1, nspins
CALL matrix_qs_to_almo(matrix_qs_ks(ispin)%matrix, &
matrix_ks(ispin), &
mat_distr_aos, .FALSE.)
CALL dbcsr_filter(matrix_ks(ispin), eps_filter)
ENDDO
CALL timestop(handle)
END SUBROUTINE almo_dm_to_almo_ks
! **************************************************************************************************
!> \brief update qs_env total energy
@ -845,9 +969,6 @@ CONTAINS
! fifth, communicate all list data
ALLOCATE (global_list(global_list_length))
!WRITE(*,*) "LENGTH: ", list_length_cpu
!WRITE(*,*) "OFFSET: ", list_offset_cpu
!WRITE(*,*) "LOCAL: ", local_list
CALL mp_allgather(local_list, global_list, &
list_length_cpu, list_offset_cpu, GroupID)
DEALLOCATE (list_length_cpu, list_offset_cpu)
@ -920,12 +1041,9 @@ CONTAINS
! O(N) loop over domain pairs
DO row = 1, nblkrows_tot
DO col = 1, current_number_neighbors(row)
!DO col = 1, nblkcols_tot
tr = .FALSE.
iblock_row = row
iblock_col = domain_neighbor_list(row, col)
!iblock_col = col
!IF (unit_nr>0) WRITE(*,*) iblock_row, iblock_col
CALL dbcsr_get_stored_coordinates(almo_scf_env%quench_t(ispin), &
iblock_row, iblock_col, hold)
@ -1178,8 +1296,9 @@ CONTAINS
contact1_radius = cp_unit_to_cp2k(contact1_radius, "angstrom")
contact2_radius = cp_unit_to_cp2k(contact2_radius, "angstrom")
!RZK-warning the procedure is faulty for molecules: the closest contacts should be found using
! the element specific radii
!RZK-warning the procedure is faulty for molecules:
! the closest contacts should be found using
! the element specific radii
! compute inner and outer cutoff radii
r0 = almo_scf_env%quencher_r0_factor*(contact1_radius+contact2_radius)
@ -1252,12 +1371,6 @@ CONTAINS
CPABORT("number of blocks is wrong")
ENDIF
! communicate local parts of the domain map
!nNodes = dbcsr_mp_numnodes(dbcsr_distribution_mp(&
! dbcsr_distribution(almo_scf_env%quench_t(ispin))))
!GroupID = dbcsr_mp_group(dbcsr_distribution_mp(&
! dbcsr_distribution(almo_scf_env%quench_t(ispin))))
! first, communicate map sizes on the other nodes
ALLOCATE (domain_entries_cpu(nNodes), offset_for_cpu(nNodes))
CALL mp_allgather(2*domain_map_local_entries, domain_entries_cpu, GroupID)
@ -1454,7 +1567,7 @@ CONTAINS
! 6. TMP_NN1=TMP_NO1.T^(tr)=(TsiginvT^(tr)FTsiginv).T^(tr)=RFR
CALL dbcsr_multiply("N", "T", 1.0_dp, tmp_no1, almo_scf_env%matrix_t(ispin), &
0.0_dp, tmp_nn1, filter_eps=almo_scf_env%eps_filter)
CALL matrix_almo_to_qs(tmp_nn1, matrix_w(ispin)%matrix, almo_scf_env)
CALL matrix_almo_to_qs(tmp_nn1, matrix_w(ispin)%matrix, almo_scf_env%mat_distr_aos)
CALL dbcsr_release(tmp_nn1)
CALL dbcsr_release(tmp_no1)

38
src/almo_scf_types.F
View file

@ -18,8 +18,8 @@ MODULE almo_scf_types
domain_submatrix_type
USE input_constants, ONLY: &
cg_dai_yuan, cg_fletcher, cg_fletcher_reeves, cg_hager_zhang, cg_hestenes_stiefel, &
cg_liu_storey, cg_polak_ribiere, cg_zero, optimizer_diis, optimizer_pcg, prec_default, &
prec_zero
cg_liu_storey, cg_polak_ribiere, cg_zero, optimizer_diis, optimizer_pcg, &
xalmo_prec_domain, xalmo_prec_full, xalmo_prec_zero
USE kinds, ONLY: dp
#include "./base/base_uses.f90"
@ -41,6 +41,14 @@ MODULE almo_scf_types
print_optimizer_options, almo_scf_env_release, &
almo_scf_history_type
! methods that add penalty terms to the energy functional
TYPE penalty_type
REAL(KIND=dp) :: occ_vol_coeff
INTEGER :: occ_vol_method
END TYPE penalty_type
! almo-based electronic structure analysis
TYPE almo_analysis_type
@ -67,17 +75,22 @@ MODULE almo_scf_types
!INTEGER :: ndiis_q -> ndiis
REAL(KIND=dp) :: eps_error, &
eps_error_early, &
lin_search_eps_error, &
lin_search_step_size_guess
lin_search_step_size_guess, &
neglect_threshold
INTEGER :: optimizer_type ! diis, pcg, etc.
INTEGER :: preconditioner, & ! preconditioner type
conjugator, & ! conjugator type
max_iter, &
max_iter_early, &
max_iter_outer_loop, &
ndiis ! diis history length
LOGICAL :: early_stopping_on = .FALSE.
END TYPE optimizer_options_type
TYPE almo_scf_history_type
@ -199,13 +212,16 @@ MODULE almo_scf_types
LOGICAL :: s_inv_done
LOGICAL :: s_sqrt_done
REAL(KIND=dp) :: almo_scf_energy
LOGICAL :: orthogonal_basis, fixed_mu
LOGICAL :: orthogonal_basis, fixed_mu
LOGICAL :: return_orthogonalized_mos
! Controls for the SCF procedure
REAL(KIND=dp) :: eps_filter
INTEGER :: xalmo_trial_wf
INTEGER :: almo_scf_guess
REAL(KIND=dp) :: eps_prev_guess
INTEGER :: order_lanczos
REAL(KIND=dp) :: matrix_iter_eps_error_factor
REAL(KIND=dp) :: eps_lanczos
INTEGER :: max_iter_lanczos
REAL(KIND=dp) :: mixing_fraction
@ -214,6 +230,8 @@ MODULE almo_scf_types
INTEGER :: almo_update_algorithm
! SCF procedure for the quenched ALMOs (xALMOs)
INTEGER :: xalmo_update_algorithm
! mo overlap inversion algorithm
INTEGER :: sigma_inv_algorithm
! ALMO SCF delocalization control
LOGICAL :: perturbative_delocalization
@ -328,14 +346,16 @@ MODULE almo_scf_types
INTEGER, DIMENSION(:), ALLOCATABLE :: cpu_of_domain
! Options for various optimizers collected neatly
! Options for various subsection options collected neatly
TYPE(almo_analysis_type) :: almo_analysis
TYPE(penalty_type) :: penalty
! Options for various optimizers collected neatly
TYPE(optimizer_options_type) :: opt_block_diag_diis
TYPE(optimizer_options_type) :: opt_block_diag_pcg
TYPE(optimizer_options_type) :: opt_xalmo_diis
TYPE(optimizer_options_type) :: opt_xalmo_pcg
TYPE(optimizer_options_type) :: opt_xalmo_newton_pcg_solver
TYPE(optimizer_options_type) :: opt_k_pcg
! keywords that control electron delocalization treatment
@ -433,10 +453,12 @@ CONTAINS
optimizer%max_iter_outer_loop
SELECT CASE (optimizer%preconditioner)
CASE (prec_zero)
CASE (xalmo_prec_zero)
prec_string = "NONE"
CASE (prec_default)
prec_string = "0.5 H + 0.5 S"
CASE (xalmo_prec_domain)
prec_string = "0.5 KS + 0.5 S, DOMAINS"
CASE (xalmo_prec_full)
prec_string = "0.5 KS + 0.5 S, FULL"
END SELECT
WRITE (unit_nr, '(T4,A,T48,A33)') "preconditioner:", TRIM(prec_string)

19
src/common/bibliography.F
View file

@ -82,7 +82,7 @@ MODULE bibliography
Bates2013, Andermatt2016, Zhu2016, Schuett2016, Lu2004, &
Becke1988b, Migliore2009, Mavros2015, Holmberg2017, Marek2014, &
Stoychev2016, Futera2017, Bailey2006, Papior2017, Lehtola2018, &
Brieuc2016, Barca2018, Huang2011, Heaton_Burgess2007
Brieuc2016, Barca2018, Scheiber2018, Huang2011, Heaton_Burgess2007
CONTAINS
@ -3863,6 +3863,23 @@ CONTAINS
"ER"), &
DOI="10.1103/PhysRevLett.98.256401")
CALL add_reference(key=Scheiber2018, ISI_record=s2a( &
"AU Scheiber, H", &
" Shi, Y", &
" Khaliullin, RZ", &
"AF Scheiber, Hayden", &
" Shi, Yifei", &
" Khaliullin, Rustam Z.", &
"TI Compact orbitals enable low-cost linear-scaling ab initio molecular dynamics for weakly-interacting systems", &
"SO The Journal of Chemical Physics", &
"PD JUN 21", &
"PY 2018", &
"VL 148", &
"AR 231103", &
"DI 10.1063/1.5029939", &
"ER"), &
DOI="10.1063/1.5029939")
END SUBROUTINE add_all_references
END MODULE bibliography

22
src/input_constants.F
View file

@ -891,7 +891,11 @@ MODULE input_constants
INTEGER, PARAMETER, PUBLIC :: almo_scf_dm_sign = 1, &
almo_scf_diag = 2, &
almo_scf_pcg = 3
almo_scf_pcg = 3, &
almo_scf_skip = 4
INTEGER, PARAMETER, PUBLIC :: almo_occ_vol_penalty_none = 0, &
almo_occ_vol_penalty_lndet = 1
! optimizer parameters
INTEGER, PARAMETER, PUBLIC :: cg_zero = 0, &
@ -903,13 +907,16 @@ MODULE input_constants
cg_dai_yuan = 6, &
cg_hager_zhang = 7
INTEGER, PARAMETER, PUBLIC :: optimizer_diis = 1, &
optimizer_pcg = 2
INTEGER, PARAMETER, PUBLIC :: prec_zero = 0, &
prec_default = -1, &
prec_ks_plus_s = 4
optimizer_pcg = 2, &
optimizer_lin_eq_pcg = 3
INTEGER, PARAMETER, PUBLIC :: xalmo_prec_zero = 0, &
xalmo_prec_domain = 1, &
xalmo_prec_full = 2
INTEGER, PARAMETER, PUBLIC :: xalmo_case_block_diag = 0, &
xalmo_case_fully_deloc = 1, &
xalmo_case_normal = -1
INTEGER, PARAMETER, PUBLIC :: xalmo_trial_simplex = 0, &
xalmo_trial_r0_out = 1
! parameters for CT methods
INTEGER, PARAMETER, PUBLIC :: tensor_orthogonal = 1, &
@ -919,6 +926,11 @@ MODULE input_constants
virt_occ_size = 3, &
virt_number = 4
! spd matrix inversion algorithm
INTEGER, PARAMETER, PUBLIC :: spd_inversion_ls_hotelling = 0, &
spd_inversion_dense_cholesky = 1, &
spd_inversion_ls_taylor = 2
! some MP2 parameters
INTEGER, PARAMETER, PUBLIC :: mp2_method_none = 0, &
mp2_method_laplace = 2, &

288
src/input_cp2k_almo.F
View file

@ -10,15 +10,20 @@
MODULE input_cp2k_almo
USE bibliography, ONLY: Khaliullin2007,&
Khaliullin2008,&
Khaliullin2013
Khaliullin2013,&
Scheiber2018
USE cp_output_handling, ONLY: cp_print_key_section_create,&
low_print_level
USE input_constants, ONLY: &
almo_deloc_none, almo_deloc_scf, almo_deloc_x, almo_deloc_x_then_scf, &
almo_deloc_xalmo_1diag, almo_deloc_xalmo_scf, almo_deloc_xalmo_x, almo_frz_crystal, &
almo_frz_none, almo_scf_diag, almo_scf_pcg, atomic_guess, cg_dai_yuan, cg_fletcher, &
cg_fletcher_reeves, cg_hager_zhang, cg_hestenes_stiefel, cg_liu_storey, cg_polak_ribiere, &
cg_zero, molecular_guess, optimizer_diis, optimizer_pcg, prec_default, prec_zero
almo_frz_none, almo_occ_vol_penalty_lndet, almo_occ_vol_penalty_none, almo_scf_diag, &
almo_scf_pcg, almo_scf_skip, atomic_guess, cg_dai_yuan, cg_fletcher, cg_fletcher_reeves, &
cg_hager_zhang, cg_hestenes_stiefel, cg_liu_storey, cg_polak_ribiere, cg_zero, &
molecular_guess, optimizer_diis, optimizer_lin_eq_pcg, optimizer_pcg, &
spd_inversion_dense_cholesky, spd_inversion_ls_hotelling, spd_inversion_ls_taylor, &
xalmo_prec_domain, xalmo_prec_full, xalmo_prec_zero, xalmo_trial_r0_out, &
xalmo_trial_simplex
USE input_keyword_types, ONLY: keyword_create,&
keyword_release,&
keyword_type
@ -39,6 +44,7 @@ MODULE input_cp2k_almo
INTEGER, PARAMETER, PRIVATE :: optimizer_block_diagonal_diis = 1
INTEGER, PARAMETER, PRIVATE :: optimizer_block_diagonal_pcg = 2
INTEGER, PARAMETER, PRIVATE :: optimizer_xalmo_pcg = 3
INTEGER, PARAMETER, PRIVATE :: optimizer_newton_pcg_solver = 5
PUBLIC :: create_almo_scf_section
@ -65,8 +71,8 @@ CONTAINS
description="Settings for a class of efficient linear scaling methods based "// &
"on absolutely localized orbitals"// &
" (ALMOs). ALMO methods are currently restricted to closed-shell molecular systems.", &
n_keywords=7, n_subsections=3, repeats=.FALSE., &
citations=(/Khaliullin2013/))
n_keywords=9, n_subsections=5, repeats=.FALSE., &
citations=(/Khaliullin2013, Scheiber2018/))
NULLIFY (keyword)
@ -89,6 +95,25 @@ CONTAINS
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="MO_OVERLAP_INV_ALG", &
description="Algorithm to invert MO overlap matrix.", &
usage="MO_OVERLAP_INV_ALG LS_HOTELLING", &
default_i_val=spd_inversion_ls_hotelling, &
enum_c_vals=s2a("LS_HOTELLING", "LS_TAYLOR", "DENSE_CHOLESKY"), &
enum_desc=s2a("Linear scaling iterative Hotelling algorithm. Fast for large sparse matrices.", &
"Linear scaling algorithm based on Taylor expansion of (A+B)^(-1).", &
"Stable but dense Cholesky algorithm. Cubically scaling."), &
enum_i_vals=(/spd_inversion_ls_hotelling, spd_inversion_ls_taylor, spd_inversion_dense_cholesky/))
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
!CALL keyword_create(keyword, name="STOP_SCF_EARLY",&
! description="Stops SCF using EPS_ERROR_EARLY or MAX_ITER_EARLY", &
! usage="STOP_SCF_EARLY .TRUE.", default_l_val=.FALSE.,&
! lone_keyword_l_val=.TRUE.)
!CALL section_add_keyword(section,keyword)
!CALL keyword_release(keyword)
CALL keyword_create(keyword, name="ALMO_EXTRAPOLATION_ORDER", &
description="Number of previous states used for the ASPC extrapolation of the ALMO "// &
"initial guess. 0 implies that the guess is given by ALMO_SCF_GUESS at each step.", &
@ -99,7 +124,7 @@ CONTAINS
CALL keyword_create(keyword, name="XALMO_EXTRAPOLATION_ORDER", &
description="Number of previous states used for the ASPC extrapolation of the initial guess "// &
"for the delocalization correction.", &
usage="XALMO_EXTRAPOLATION_ORDER 1", default_i_val=3)
usage="XALMO_EXTRAPOLATION_ORDER 1", default_i_val=0)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
@ -108,14 +133,28 @@ CONTAINS
usage="ALMO_ALGORITHM DIAG", &
default_i_val=almo_scf_diag, &
!enum_c_vals=s2a("DIAG", "DM_SIGN","PCG"),&
enum_c_vals=s2a("DIAG", "PCG"), &
!enum_c_vals=s2a("DIAG", "PCG"), &
enum_c_vals=s2a("DIAG", "PCG", "SKIP"), &
enum_desc=s2a("DIIS-accelerated diagonalization controlled by ALMO_OPTIMIZER_DIIS. "// &
"Recommended for large systems containing small fragments.", &
!"Update the density matrix using linear scaling routines. "//&
!"Recommended if large fragments are present.",&
"Energy minimization with a PCG algorithm controlled by ALMO_OPTIMIZER_PCG."), &
"Energy minimization with a PCG algorithm controlled by ALMO_OPTIMIZER_PCG.", &
"Skip optimization of block-diagonal ALMOs."), &
!enum_i_vals=(/almo_scf_diag,almo_scf_dm_sign,almo_scf_pcg/),&
enum_i_vals=(/almo_scf_diag, almo_scf_pcg/))
enum_i_vals=(/almo_scf_diag, almo_scf_pcg, almo_scf_skip/))
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="XALMO_TRIAL_WF", &
description="Determines the form of the trial XALMOs.", &
usage="XALMO_TRIAL_WF SIMPLE", &
default_i_val=xalmo_trial_r0_out, &
enum_c_vals=s2a("SIMPLE", "PROJECT_R0_OUT"), &
enum_desc=s2a("Straightforward AO-basis expansion.", &
"Block-diagonal ALMOs plus the XALMO term projected onto the unoccupied "// &
"ALMO-subspace."), &
enum_i_vals=(/xalmo_trial_simplex, xalmo_trial_r0_out/))
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
@ -150,6 +189,14 @@ CONTAINS
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="RETURN_ORTHOGONALIZED_MOS", &
description="Orthogonalize final ALMOs before they are returned"// &
"to Quickstep (i.e. for calculation of properties)", &
usage="RETURN_ORTHOGONALIZED_MOS .TRUE.", default_l_val=.TRUE., &
lone_keyword_l_val=.TRUE.)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
!CALL keyword_create(keyword, name="DELOCALIZE_EPS_ITER",&
! description="Obsolete and to be deleted: use EPS_ERROR in XALMO_OPTIMIZER_PCG",&
! usage="DELOCALIZE_EPS_ITER 1.e-5", default_r_val=1.e-5_dp)
@ -542,6 +589,16 @@ CONTAINS
CALL section_add_subsection(section, subsection)
CALL section_release(subsection)
NULLIFY (subsection)
CALL create_penalty_section(subsection)
CALL section_add_subsection(section, subsection)
CALL section_release(subsection)
NULLIFY (subsection)
CALL create_matrix_iterate_section(subsection)
CALL section_add_subsection(section, subsection)
CALL section_release(subsection)
NULLIFY (subsection)
CALL create_almo_analysis_section(subsection)
CALL section_add_subsection(section, subsection)
@ -559,7 +616,7 @@ CONTAINS
!> 2014.10 fully integrated [Rustam Z Khaliullin]
!> \author Rustam Z Khaliullin
! **************************************************************************************************
SUBROUTINE create_optimizer_section(section, optimizer_id)
RECURSIVE SUBROUTINE create_optimizer_section(section, optimizer_id)
TYPE(section_type), POINTER :: section
INTEGER, INTENT(IN) :: optimizer_id
@ -569,6 +626,7 @@ CONTAINS
INTEGER :: optimizer_type
TYPE(keyword_type), POINTER :: keyword
TYPE(section_type), POINTER :: subsection
CPASSERT(.NOT. ASSOCIATED(section))
NULLIFY (section)
@ -578,18 +636,27 @@ CONTAINS
CASE (optimizer_block_diagonal_diis)
CALL section_create(section, "ALMO_OPTIMIZER_DIIS", &
description="Controls the iterative DIIS-accelerated optimization of block-diagonal ALMOs.", &
n_keywords=3, n_subsections=0, repeats=.FALSE.)
n_keywords=5, n_subsections=0, repeats=.FALSE.)
optimizer_type = optimizer_diis
CASE (optimizer_block_diagonal_pcg)
CALL section_create(section, "ALMO_OPTIMIZER_PCG", &
description="Controls the PCG optimization of block-diagonal ALMOs.", &
n_keywords=6, n_subsections=0, repeats=.FALSE.)
n_keywords=9, n_subsections=0, repeats=.FALSE.)
optimizer_type = optimizer_pcg
CASE (optimizer_xalmo_pcg)
CALL section_create(section, "XALMO_OPTIMIZER_PCG", &
description="Controls the PCG optimization of extended ALMOs.", &
n_keywords=6, n_subsections=0, repeats=.FALSE.)
n_keywords=10, n_subsections=1, repeats=.FALSE.)
NULLIFY (subsection)
CALL create_optimizer_section(subsection, optimizer_newton_pcg_solver)
CALL section_add_subsection(section, subsection)
CALL section_release(subsection)
optimizer_type = optimizer_pcg
CASE (optimizer_newton_pcg_solver)
CALL section_create(section, "XALMO_NEWTON_PCG_SOLVER", &
description="Controls an iterative solver of the Newton-Raphson linear equation.", &
n_keywords=4, n_subsections=0, repeats=.FALSE.)
optimizer_type = optimizer_lin_eq_pcg
CASE DEFAULT
CPABORT("No default values allowed")
END SELECT
@ -609,6 +676,21 @@ CONTAINS
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
! add common keywords
CALL keyword_create(keyword, name="MAX_ITER_EARLY", &
description="Maximum number of iterations for truncated SCF "// &
"(e.g. Langevin-corrected MD). Negative values mean that this keyword is not used.", &
usage="MAX_ITER_EARLY 5", default_i_val=-1)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="EPS_ERROR_EARLY", &
description="Target value of the MAX norm of the error for truncated SCF "// &
"(e.g. Langevin-corrected MD). Negative values mean that this keyword is not used.", &
usage="EPS_ERROR_EARLY 1.E-2", default_r_val=-1.0_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
! add keywords specific to each type
SELECT CASE (optimizer_type)
CASE (optimizer_diis)
@ -620,19 +702,50 @@ CONTAINS
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CASE (optimizer_pcg)
CASE (optimizer_pcg, optimizer_lin_eq_pcg)
CALL keyword_create(keyword, name="LIN_SEARCH_EPS_ERROR", &
description="Target value of the gradient norm during the linear search", &
usage="LIN_SEARCH_EPS_ERROR 1.E-2", default_r_val=1.0E-3_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
IF (optimizer_type .EQ. optimizer_pcg) THEN
CALL keyword_create(keyword, name="LIN_SEARCH_STEP_SIZE_GUESS", &
description="The size of the first step in the linear search", &
usage="LIN_SEARCH_STEP_SIZE_GUESS 0.1", default_r_val=1.0_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="LIN_SEARCH_EPS_ERROR", &
description="Target value of the gradient norm during the linear search", &
usage="LIN_SEARCH_EPS_ERROR 1.E-2", default_r_val=1.0E-3_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="LIN_SEARCH_STEP_SIZE_GUESS", &
description="The size of the first step in the linear search", &
usage="LIN_SEARCH_STEP_SIZE_GUESS 0.1", default_r_val=1.0_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="CONJUGATOR", &
description="Various methods to compute step directions in the PCG optimization", &
usage="CONJUGATOR POLAK_RIBIERE", &
default_i_val=cg_hager_zhang, &
enum_c_vals=s2a("ZERO", "POLAK_RIBIERE", "FLETCHER_REEVES", &
"HESTENES_STIEFEL", "FLETCHER", "LIU_STOREY", "DAI_YUAN", "HAGER_ZHANG"), &
enum_desc=s2a("Steepest descent", "Polak and Ribiere", &
"Fletcher and Reeves", "Hestenes and Stiefel", &
"Fletcher (Conjugate descent)", "Liu and Storey", &
"Dai and Yuan", "Hager and Zhang"), &
enum_i_vals=(/cg_zero, cg_polak_ribiere, cg_fletcher_reeves, &
cg_hestenes_stiefel, cg_fletcher, cg_liu_storey, &
cg_dai_yuan, cg_hager_zhang/))
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="PRECOND_FILTER_THRESHOLD", &
description="Select eigenvalues of the preconditioner "// &
"that are smaller than the threshold and project out the "// &
"corresponding eigenvectors from the gradient. No matter "// &
"how large the threshold is the maximum number of projected "// &
"eienvectors for a fragment equals to the number of occupied "// &
"orbitals of fragment's neighbors.", &
usage="PRECOND_FILTER_THRESHOLD 0.1", default_r_val=0.5_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
ENDIF
CALL keyword_create(keyword, name="MAX_ITER_OUTER_LOOP", &
description="Maximum number of iterations in the outer loop. "// &
@ -644,27 +757,16 @@ CONTAINS
CALL keyword_create(keyword, name="PRECONDITIONER", &
description="Select a preconditioner for the conjugate gradient optimization", &
usage="PRECONDITIONER NONE", &
default_i_val=-1, &
enum_c_vals=s2a("DEFAULT", "NONE"), &
enum_desc=s2a("Default preconditioner", "Do not use preconditioner"), &
enum_i_vals=(/prec_default, prec_zero/))
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="CONJUGATOR", &
description="Various methods to compute step directions in the PCG optimization", &
usage="CONJUGATOR POLAK_RIBIERE", &
default_i_val=cg_hager_zhang, &
enum_c_vals=s2a("ZERO", "POLAK_RIBIERE", "FLETCHER_REEVES", &
"HESTENES_STIEFEL", "FLETCHER", "LIU_STOREY", "DAI_YUAN", "HAGER_ZHANG"), &
enum_desc=s2a("Steepest descent", "Polak and Ribiere", &
"Fletcher and Reeves", "Hestenes and Stiefel", &
"Fletcher (Conjugate descent)", "Liu and Storey", &
"Dai and Yuan", "Hager and Zhang"), &
enum_i_vals=(/cg_zero, cg_polak_ribiere, cg_fletcher_reeves, &
cg_hestenes_stiefel, cg_fletcher, cg_liu_storey, &
cg_dai_yuan, cg_hager_zhang/))
usage="PRECONDITIONER DOMAIN", &
default_i_val=1, &
enum_c_vals=s2a("NONE", "DEFAULT", "DOMAIN", "FULL"), &
enum_desc=s2a("Do not use preconditioner", &
"Same as DOMAIN preconditioner", &
"Invert preconditioner domain-by-domain."// &
" The main component of the linear scaling algorithm", &
"Solve linear equations step=-H.grad on the entire space"), &
enum_i_vals=(/xalmo_prec_zero, xalmo_prec_domain, &
xalmo_prec_domain, xalmo_prec_full/))
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
@ -672,6 +774,98 @@ CONTAINS
END SUBROUTINE create_optimizer_section
! **************************************************************************************************
!> \brief The section controls iterative matrix operations like SQRT or inverse
!> \param section ...
!> \par History
!> 2017.05 created [Rustam Z Khaliullin]
!> \author Rustam Z Khaliullin
! **************************************************************************************************
SUBROUTINE create_matrix_iterate_section(section)
TYPE(section_type), POINTER :: section
TYPE(keyword_type), POINTER :: keyword
CPASSERT(.NOT. ASSOCIATED(section))
NULLIFY (section)
CALL section_create(section, "MATRIX_ITERATE", &
description="Controls linear scaling iterative procedure on matrices: inversion, sqrti, etc. "// &
"High-order Lanczos accelerates convergence provided it can estimate the eigenspectrum correctly.", &
n_keywords=4, n_subsections=0, repeats=.FALSE.)
NULLIFY (keyword)
CALL keyword_create(keyword, name="EPS_TARGET_FACTOR", &
description="Multiplication factor that determines acceptable error in the iterative procedure. "// &
"Acceptable error = EPS_TARGET_FACTOR * EPS_FILTER", &
usage="EPS_TARGET_FACTOR 100.0", default_r_val=10.0_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="EPS_LANCZOS", &
description="Threshold for Lanczos eigenvalue estimation.", &
usage="EPS_LANCZOS 1.0E-4", default_r_val=1.0E-3_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="ORDER_LANCZOS", &
description="Order of the Lanczos estimator. Use 0 to turn off. Do not use 1.", &
usage="ORDER_LANCZOS 5", default_i_val=3)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="MAX_ITER_LANCZOS", &
description="Maximum number of Lanczos iterations.", &
usage="MAX_ITER_LANCZOS 64", default_i_val=128)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
END SUBROUTINE create_matrix_iterate_section
! **************************************************************************************************
!> \brief The section controls penalty methods
!> \param section ...
!> \par History
!> 2018.01 created [Rustam Z Khaliullin]
!> \author Rustam Z Khaliullin
! **************************************************************************************************
SUBROUTINE create_penalty_section(section)
TYPE(section_type), POINTER :: section
TYPE(keyword_type), POINTER :: keyword
CPASSERT(.NOT. ASSOCIATED(section))
NULLIFY (section)
CALL section_create(section, "PENALTY", &
description="Add penalty terms to the energy functional.", &
n_keywords=2, n_subsections=0, repeats=.FALSE.)
NULLIFY (keyword)
CALL keyword_create( &
keyword, name="OCCUPIED_VOLUME_PENALTY_METHOD", &
description="Penalty that prevents nonorthogonal orbitals from becoming linear dependent.", &
usage="OCCUPIED_VOLUME_PENALTY_METHOD LNDET", &
default_i_val=almo_occ_vol_penalty_none, &
enum_c_vals=s2a("NONE", "LNDET"), &
enum_desc=s2a("Do not use penalties", &
"Use -coeff*ln(det(MO-overlap)) term."), &
enum_i_vals=(/almo_occ_vol_penalty_none, almo_occ_vol_penalty_lndet/))
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="OCCUPIED_VOLUME_PENALTY_COEFF", &
description="Multiplication factor that determines the strength of the penalty term.", &
usage="OCCUPIED_VOLUME_PENALTY_COEFF 10.0", default_r_val=1.0_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
END SUBROUTINE create_penalty_section
! **************************************************************************************************
!> \brief The section controls electronic structure analysis based on ALMOs
!> \param section ...
@ -694,7 +888,7 @@ CONTAINS
CALL section_create(section, "ANALYSIS", &
description="Controls electronic structure analysis based on ALMOs and XALMOs.", &
n_keywords=2, n_subsections=0, repeats=.FALSE., &
n_keywords=2, n_subsections=1, repeats=.FALSE., &
citations=(/Khaliullin2007, Khaliullin2008/))
NULLIFY (keyword)

594
src/iterate_matrix.F
View file

@ -10,20 +10,16 @@
! **************************************************************************************************
MODULE iterate_matrix
USE arnoldi_api, ONLY: arnoldi_data_type,&
arnoldi_ev,&
arnoldi_extremal,&
deallocate_arnoldi_data,&
get_selected_ritz_val,&
setup_arnoldi_data
arnoldi_extremal
USE cp_log_handling, ONLY: cp_get_default_logger,&
cp_logger_get_default_unit_nr,&
cp_logger_type
USE dbcsr_api, ONLY: &
dbcsr_add, dbcsr_add_on_diag, dbcsr_copy, dbcsr_create, dbcsr_filter, &
dbcsr_frobenius_norm, dbcsr_gershgorin_norm, dbcsr_get_info, dbcsr_get_matrix_type, &
dbcsr_get_occupation, dbcsr_multiply, dbcsr_norm, dbcsr_norm_maxabsnorm, dbcsr_p_type, &
dbcsr_release, dbcsr_scale, dbcsr_set, dbcsr_trace, dbcsr_transposed, dbcsr_type, &
dbcsr_type_no_symmetry
dbcsr_add, dbcsr_add_on_diag, dbcsr_copy, dbcsr_create, dbcsr_desymmetrize, dbcsr_filter, &
dbcsr_frobenius_norm, dbcsr_gershgorin_norm, dbcsr_get_diag, dbcsr_get_info, &
dbcsr_get_matrix_type, dbcsr_get_occupation, dbcsr_multiply, dbcsr_norm, &
dbcsr_norm_maxabsnorm, dbcsr_p_type, dbcsr_release, dbcsr_scale, dbcsr_set, &
dbcsr_set_diag, dbcsr_trace, dbcsr_transposed, dbcsr_type, dbcsr_type_no_symmetry
USE kinds, ONLY: dp,&
int_8
USE machine, ONLY: m_flush,&
@ -43,51 +39,246 @@ MODULE iterate_matrix
END INTERFACE
PUBLIC :: invert_Hotelling, matrix_sign_Newton_Schulz, matrix_sqrt_Newton_Schulz, &
purify_mcweeny
purify_mcweeny, invert_Taylor, determinant
CONTAINS
! *****************************************************************************
!> \brief Computes the determinant of a symmetric positive definite matrix
!> using the trace of the matrix logarithm via Mercator series:
!> det(A) = det(S)det(I+X)det(S), where S=diag(sqrt(Aii),..,sqrt(Ann))
!> det(I+X) = Exp(Trace(Ln(I+X)))
!> Ln(I+X) = X - X^2/2 + X^3/3 - X^4/4 + ..
!> The series converges only if the Frobenius norm of X is less than 1.
!> If it is more than one we compute (recursevily) the determinant of
!> the square root of (I+X).
!> \param matrix ...
!> \param det - determinant
!> \param threshold ...
!> \par History
!> 2015.04 created [Rustam Z Khaliullin]
!> \author Rustam Z. Khaliullin
! **************************************************************************************************
!> \brief invert a symmetric positive definite matrix by Hotelling's method
!> explicit symmetrization makes this code not suitable for other matrix types
!> Currently a bit messy with the options, to to be cleaned soon
RECURSIVE SUBROUTINE determinant(matrix, det, threshold)
TYPE(dbcsr_type), INTENT(INOUT) :: matrix
REAL(KIND=dp), INTENT(INOUT) :: det
REAL(KIND=dp), INTENT(IN) :: threshold
CHARACTER(LEN=*), PARAMETER :: routineN = 'determinant', routineP = moduleN//':'//routineN
INTEGER :: handle, i, max_iter_lanczos, nsize, &
order_lanczos, sign_iter, unit_nr
INTEGER(KIND=int_8) :: flop1, flop2
INTEGER, SAVE :: recursion_depth = 0
REAL(KIND=dp) :: det0, eps_lanczos, frobnorm, maxnorm, &
occ_matrix, t1, t2, trace
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: diagonal
TYPE(cp_logger_type), POINTER :: logger
TYPE(dbcsr_type) :: tmp1, tmp2, tmp3
CALL timeset(routineN, handle)
! get a useful output_unit
logger => cp_get_default_logger()
IF (logger%para_env%mepos == logger%para_env%source) THEN
unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.)
ELSE
unit_nr = -1
ENDIF
! Note: tmp1 and tmp2 have the same matrix type as the
! initial matrix (tmp3 does not have symmetry constraints)
! this might lead to uninteded results with anti-symmetric
! matrices
CALL dbcsr_create(tmp1, template=matrix, &
matrix_type=dbcsr_type_no_symmetry)
CALL dbcsr_create(tmp2, template=matrix, &
matrix_type=dbcsr_type_no_symmetry)
CALL dbcsr_create(tmp3, template=matrix, &
matrix_type=dbcsr_type_no_symmetry)
! compute the product of the diagonal elements
CALL dbcsr_get_info(matrix, nfullrows_total=nsize)
ALLOCATE (diagonal(nsize))
CALL dbcsr_get_diag(matrix, diagonal)
det = PRODUCT(diagonal)
! create diagonal SQRTI matrix
diagonal(:) = 1.0_dp/(SQRT(diagonal(:)))
!ROLL CALL dbcsr_copy(tmp1,matrix)
CALL dbcsr_desymmetrize(matrix, tmp1)
CALL dbcsr_set(tmp1, 0.0_dp)
CALL dbcsr_set_diag(tmp1, diagonal)
CALL dbcsr_filter(tmp1, threshold)
DEALLOCATE (diagonal)
! normalize the main diagonal, off-diagonal elements are scaled to
! make the norm of the matrix less than 1
CALL dbcsr_multiply("N", "N", 1.0_dp, &
matrix, &
tmp1, &
0.0_dp, tmp3, &
filter_eps=threshold)
CALL dbcsr_multiply("N", "N", 1.0_dp, &
tmp1, &
tmp3, &
0.0_dp, tmp2, &
filter_eps=threshold)
! subtract the main diagonal to create matrix X
CALL dbcsr_add_on_diag(tmp2, -1.0_dp)
frobnorm = dbcsr_frobenius_norm(tmp2)
IF (unit_nr > 0) THEN
IF (recursion_depth .EQ. 0) THEN
WRITE (unit_nr, '()')
ELSE
WRITE (unit_nr, '(T6,A28,1X,I15)') &
"Recursive iteration:", recursion_depth
ENDIF
WRITE (unit_nr, '(T6,A28,1X,F15.10)') &
"Frobenius norm:", frobnorm
CALL m_flush(unit_nr)
ENDIF
IF (frobnorm .GE. 1.0_dp) THEN
CALL dbcsr_add_on_diag(tmp2, 1.0_dp)
! these controls should be provided as input
order_lanczos = 3
eps_lanczos = 1.0E-4_dp
max_iter_lanczos = 40
CALL matrix_sqrt_Newton_Schulz( &
tmp3, & ! output sqrt
tmp1, & ! output sqrti
tmp2, & ! input original
threshold=threshold, &
order=order_lanczos, &
eps_lanczos=eps_lanczos, &
max_iter_lanczos=max_iter_lanczos)
recursion_depth = recursion_depth+1
CALL determinant(tmp3, det0, threshold)
recursion_depth = recursion_depth-1
det = det*det0*det0
ELSE
! create accumulator
CALL dbcsr_copy(tmp1, tmp2)
! re-create to make use of symmetry
!ROLL CALL dbcsr_create(tmp3,template=matrix)
IF (unit_nr > 0) WRITE (unit_nr, *)
! initialize the sign of the term
sign_iter = -1
DO i = 1, 100
t1 = m_walltime()
! multiply X^i by X
! note that the first iteration evaluates X^2
! because the trace of X^1 is zero by construction
CALL dbcsr_multiply("N", "N", 1.0_dp, tmp1, tmp2, &
0.0_dp, tmp3, &
filter_eps=threshold, &
flop=flop1)
CALL dbcsr_copy(tmp1, tmp3)
! get trace
CALL dbcsr_trace(tmp1, trace)
trace = trace*sign_iter/(1.0_dp*(i+1))
sign_iter = -sign_iter
! update the determinant
det = det*EXP(trace)
occ_matrix = dbcsr_get_occupation(tmp1)
CALL dbcsr_norm(tmp1, &
dbcsr_norm_maxabsnorm, norm_scalar=maxnorm)
t2 = m_walltime()
IF (unit_nr > 0) THEN
WRITE (unit_nr, '(T6,A,1X,I3,1X,F7.5,F16.10,F10.3,F11.3)') &
"Determinant iter", i, occ_matrix, &
det, t2-t1, &
(flop1+flop2)/(1.0E6_dp*MAX(0.001_dp, t2-t1))
CALL m_flush(unit_nr)
ENDIF
! exit if the trace is close to zero
IF (maxnorm < threshold) EXIT
ENDDO ! end iterations
IF (unit_nr > 0) THEN
WRITE (unit_nr, '()')
CALL m_flush(unit_nr)
ENDIF
ENDIF ! decide to do sqrt or not
IF (unit_nr > 0) THEN
IF (recursion_depth .EQ. 0) THEN
WRITE (unit_nr, '(T6,A28,1X,F15.10)') &
"Final determinant:", det
WRITE (unit_nr, '()')
ELSE
WRITE (unit_nr, '(T6,A28,1X,F15.10)') &
"Recursive determinant:", det
ENDIF
CALL m_flush(unit_nr)
ENDIF
CALL dbcsr_release(tmp1)
CALL dbcsr_release(tmp2)
CALL dbcsr_release(tmp3)
CALL timestop(handle)
END SUBROUTINE determinant
! **************************************************************************************************
!> \brief invert a symmetric positive definite diagonally dominant matrix
!> \param matrix_inverse ...
!> \param matrix ...
!> \param threshold convergence threshold nased on the max abs
!> \param use_inv_as_guess logical whether input can be used as guess for inverse
!> \param norm_convergence convergence threshold for the 2-norm, useful for approximate solutions
!> \param filter_eps filter_eps for matrix multiplications, if not passed nothing is filteres
!> \param accelerator_order ...
!> \param max_iter_lanczos ...
!> \param eps_lanczos ...
!> \param silent ...
!> \par History
!> 2010.10 created [Joost VandeVondele]
!> 2011.10 guess option added [Rustam Z Khaliullin]
!> \author Joost VandeVondele
! **************************************************************************************************
SUBROUTINE invert_Hotelling(matrix_inverse, matrix, threshold, use_inv_as_guess, &
norm_convergence, filter_eps, silent)
SUBROUTINE invert_Taylor(matrix_inverse, matrix, threshold, use_inv_as_guess, &
norm_convergence, filter_eps, accelerator_order, &
max_iter_lanczos, eps_lanczos, silent)
TYPE(dbcsr_type), INTENT(INOUT), TARGET :: matrix_inverse, matrix
REAL(KIND=dp), INTENT(IN) :: threshold
LOGICAL, INTENT(IN), OPTIONAL :: use_inv_as_guess
REAL(KIND=dp), INTENT(IN), OPTIONAL :: norm_convergence, filter_eps
INTEGER, INTENT(IN), OPTIONAL :: accelerator_order, max_iter_lanczos
REAL(KIND=dp), INTENT(IN), OPTIONAL :: eps_lanczos
LOGICAL, INTENT(IN), OPTIONAL :: silent
CHARACTER(LEN=*), PARAMETER :: routineN = 'invert_Hotelling', &
routineP = moduleN//':'//routineN
CHARACTER(LEN=*), PARAMETER :: routineN = 'invert_Taylor', routineP = moduleN//':'//routineN
INTEGER :: handle, i, nrow, unit_nr
INTEGER(KIND=int_8) :: flop1, flop2
LOGICAL :: use_inv_guess
REAL(KIND=dp) :: convergence, frob_matrix, &
gershgorin_norm, max_ev, &
maxnorm_matrix, min_ev, occ_matrix, &
t1, t2
TYPE(arnoldi_data_type) :: my_arnoldi
INTEGER :: accelerator_type, handle, i, &
my_max_iter_lanczos, nrows, unit_nr
INTEGER(KIND=int_8) :: flop2
LOGICAL :: converged, use_inv_guess
REAL(KIND=dp) :: coeff, convergence, maxnorm_matrix, &
my_eps_lanczos, occ_matrix, t1, t2
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: p_diagonal
TYPE(cp_logger_type), POINTER :: logger
TYPE(dbcsr_p_type), DIMENSION(1) :: mymat
TYPE(dbcsr_type), TARGET :: tmp1, tmp2
! turn this off for the time being
TYPE(dbcsr_type), TARGET :: tmp1, tmp2, tmp3_sym
CALL timeset(routineN, handle)
@ -104,73 +295,301 @@ CONTAINS
convergence = threshold
IF (PRESENT(norm_convergence)) convergence = norm_convergence
accelerator_type = 0
IF (PRESENT(accelerator_order)) accelerator_type = accelerator_order
IF (accelerator_type .GT. 1) accelerator_type = 1
use_inv_guess = .FALSE.
IF (PRESENT(use_inv_as_guess)) use_inv_guess = use_inv_as_guess
my_max_iter_lanczos = 64
my_eps_lanczos = 1.0E-3_dp
IF (PRESENT(max_iter_lanczos)) my_max_iter_lanczos = max_iter_lanczos
IF (PRESENT(eps_lanczos)) my_eps_lanczos = eps_lanczos
CALL dbcsr_create(tmp1, template=matrix_inverse, matrix_type=dbcsr_type_no_symmetry)
CALL dbcsr_create(tmp2, template=matrix_inverse, matrix_type=dbcsr_type_no_symmetry)
CALL dbcsr_create(tmp3_sym, template=matrix_inverse)
CALL dbcsr_get_info(matrix, nfullrows_total=nrows)
ALLOCATE (p_diagonal(nrows))
! generate the initial guess
IF (.NOT. use_inv_guess) THEN
! initialize matrix to unity and use arnoldi to scale it into the convergence range
gershgorin_norm = dbcsr_gershgorin_norm(matrix)
frob_matrix = dbcsr_frobenius_norm(matrix)
CALL dbcsr_set(matrix_inverse, 0.0_dp)
CALL dbcsr_add_on_diag(matrix_inverse, 1.0_dp)
! everything commutes, therefor our all products will be symmetric
CALL dbcsr_create(tmp1, template=matrix_inverse)
SELECT CASE (accelerator_type)
CASE (0)
! use tmp1 to hold off-diagonal elements
CALL dbcsr_desymmetrize(matrix, tmp1)
p_diagonal(:) = 0.0_dp
CALL dbcsr_set_diag(tmp1, p_diagonal)
!CALL dbcsr_print(tmp1)
! invert the main diagonal
CALL dbcsr_get_diag(matrix, p_diagonal)
p_diagonal(:) = 1.0_dp/p_diagonal(:)
CALL dbcsr_set(matrix_inverse, 0.0_dp)
CALL dbcsr_add_on_diag(matrix_inverse, 1.0_dp)
CALL dbcsr_set_diag(matrix_inverse, p_diagonal)
CASE DEFAULT
CPABORT("Illegal accelerator order")
END SELECT
ELSE
CPABORT("Guess is NYI")
ENDIF
CALL dbcsr_multiply("N", "N", 1.0_dp, tmp1, matrix_inverse, &
0.0_dp, tmp2, filter_eps=filter_eps)
IF (unit_nr > 0) WRITE (unit_nr, *)
! scale the approximate inverse to be within the convergence radius
t1 = m_walltime()
! done with the initial guess, start iterations
converged = .FALSE.
CALL dbcsr_desymmetrize(matrix_inverse, tmp1)
coeff = 1.0_dp
DO i = 1, 100
! coeff = +/- 1
coeff = -1.0_dp*coeff
CALL dbcsr_multiply("N", "N", 1.0_dp, tmp1, tmp2, 0.0_dp, &
tmp3_sym, &
flop=flop2, filter_eps=filter_eps)
!flop=flop2)
CALL dbcsr_add(matrix_inverse, tmp3_sym, 1.0_dp, coeff)
CALL dbcsr_release(tmp1)
CALL dbcsr_create(tmp1, template=matrix_inverse, matrix_type=dbcsr_type_no_symmetry)
CALL dbcsr_desymmetrize(tmp3_sym, tmp1)
! for the convergence check
CALL dbcsr_norm(tmp3_sym, &
dbcsr_norm_maxabsnorm, norm_scalar=maxnorm_matrix)
t2 = m_walltime()
occ_matrix = dbcsr_get_occupation(matrix_inverse)
IF (unit_nr > 0) THEN
WRITE (unit_nr, '(T6,A,1X,I3,1X,F10.8,E12.3,F12.3,F13.3)') "Taylor iter", i, occ_matrix, &
maxnorm_matrix, t2-t1, &
flop2/(1.0E6_dp*MAX(0.001_dp, t2-t1))
CALL m_flush(unit_nr)
ENDIF
IF (maxnorm_matrix < convergence) THEN
converged = .TRUE.
EXIT
ENDIF
t1 = m_walltime()
ENDDO
!last convergence check
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix, matrix_inverse, 0.0_dp, tmp1, &
filter_eps=filter_eps)
CALL dbcsr_add_on_diag(tmp1, -1.0_dp)
!frob_matrix = dbcsr_frobenius_norm(tmp1)
CALL dbcsr_norm(tmp1, dbcsr_norm_maxabsnorm, norm_scalar=maxnorm_matrix)
IF (unit_nr > 0) THEN
WRITE (unit_nr, '(T6,A,E12.5)') "Final Taylor error", maxnorm_matrix
WRITE (unit_nr, '()')
CALL m_flush(unit_nr)
ENDIF
IF (maxnorm_matrix > convergence) THEN
converged = .FALSE.
IF (unit_nr > 0) THEN
WRITE (*, *) 'Final convergence check failed'
ENDIF
ENDIF
IF (.NOT. converged) THEN
CPABORT("Taylor inversion did not converge")
ENDIF
CALL dbcsr_release(tmp1)
CALL dbcsr_release(tmp2)
CALL dbcsr_release(tmp3_sym)
DEALLOCATE (p_diagonal)
CALL timestop(handle)
END SUBROUTINE invert_Taylor
! **************************************************************************************************
!> \brief invert a symmetric positive definite matrix by Hotelling's method
!> explicit symmetrization makes this code not suitable for other matrix types
!> Currently a bit messy with the options, to to be cleaned soon
!> \param matrix_inverse ...
!> \param matrix ...
!> \param threshold convergence threshold nased on the max abs
!> \param use_inv_as_guess logical whether input can be used as guess for inverse
!> \param norm_convergence convergence threshold for the 2-norm, useful for approximate solutions
!> \param filter_eps filter_eps for matrix multiplications, if not passed nothing is filteres
!> \param accelerator_order ...
!> \param max_iter_lanczos ...
!> \param eps_lanczos ...
!> \param silent ...
!> \par History
!> 2010.10 created [Joost VandeVondele]
!> 2011.10 guess option added [Rustam Z Khaliullin]
!> \author Joost VandeVondele
! **************************************************************************************************
SUBROUTINE invert_Hotelling(matrix_inverse, matrix, threshold, use_inv_as_guess, &
norm_convergence, filter_eps, accelerator_order, &
max_iter_lanczos, eps_lanczos, silent)
TYPE(dbcsr_type), INTENT(INOUT), TARGET :: matrix_inverse, matrix
REAL(KIND=dp), INTENT(IN) :: threshold
LOGICAL, INTENT(IN), OPTIONAL :: use_inv_as_guess
REAL(KIND=dp), INTENT(IN), OPTIONAL :: norm_convergence, filter_eps
INTEGER, INTENT(IN), OPTIONAL :: accelerator_order, max_iter_lanczos
REAL(KIND=dp), INTENT(IN), OPTIONAL :: eps_lanczos
LOGICAL, INTENT(IN), OPTIONAL :: silent
CHARACTER(LEN=*), PARAMETER :: routineN = 'invert_Hotelling', &
routineP = moduleN//':'//routineN
INTEGER :: accelerator_type, handle, i, &
my_max_iter_lanczos, unit_nr
INTEGER(KIND=int_8) :: flop1, flop2
LOGICAL :: arnoldi_converged, converged, &
use_inv_guess
REAL(KIND=dp) :: convergence, frob_matrix, gershgorin_norm, max_ev, maxnorm_matrix, min_ev, &
my_eps_lanczos, my_filter_eps, occ_matrix, scalingf, t1, t2
TYPE(cp_logger_type), POINTER :: logger
TYPE(dbcsr_type), TARGET :: tmp1, tmp2
!TYPE(arnoldi_data_type) :: my_arnoldi
!TYPE(dbcsr_p_type), DIMENSION(1) :: mymat
CALL timeset(routineN, handle)
logger => cp_get_default_logger()
IF (logger%para_env%mepos == logger%para_env%source) THEN
unit_nr = cp_logger_get_default_unit_nr(logger, local=.TRUE.)
ELSE
unit_nr = -1
ENDIF
IF (PRESENT(silent)) THEN
IF (silent) unit_nr = -1
END IF
convergence = threshold
IF (PRESENT(norm_convergence)) convergence = norm_convergence
accelerator_type = 1
IF (PRESENT(accelerator_order)) accelerator_type = accelerator_order
IF (accelerator_type .GT. 1) accelerator_type = 1
use_inv_guess = .FALSE.
IF (PRESENT(use_inv_as_guess)) use_inv_guess = use_inv_as_guess
my_max_iter_lanczos = 64
my_eps_lanczos = 1.0E-3_dp
IF (PRESENT(max_iter_lanczos)) my_max_iter_lanczos = max_iter_lanczos
IF (PRESENT(eps_lanczos)) my_eps_lanczos = eps_lanczos
my_filter_eps = threshold
IF (PRESENT(filter_eps)) my_filter_eps = filter_eps
! generate the initial guess
IF (.NOT. use_inv_guess) THEN
SELECT CASE (accelerator_type)
CASE (0)
gershgorin_norm = dbcsr_gershgorin_norm(matrix)
frob_matrix = dbcsr_frobenius_norm(matrix)
CALL dbcsr_set(matrix_inverse, 0.0_dp)
CALL dbcsr_add_on_diag(matrix_inverse, 1/MIN(gershgorin_norm, frob_matrix))
CASE (1)
! initialize matrix to unity and use arnoldi (below) to scale it into the convergence range
CALL dbcsr_set(matrix_inverse, 0.0_dp)
CALL dbcsr_add_on_diag(matrix_inverse, 1.0_dp)
CASE DEFAULT
CPABORT("Illegal accelerator order")
END SELECT
! everything commutes, therefore our all products will be symmetric
CALL dbcsr_create(tmp1, template=matrix_inverse)
ELSE
! It is unlikely that our guess will commute with the matrix, therefore the first product will
! be non symmetric
CALL dbcsr_create(tmp1, template=matrix_inverse, matrix_type=dbcsr_type_no_symmetry)
ENDIF
CALL dbcsr_get_info(matrix, nfullrows_total=nrow)
CALL dbcsr_create(tmp2, template=matrix_inverse)
IF (unit_nr > 0) WRITE (unit_nr, *)
! scale the approximate inverse to be within the convergence radius
t1 = m_walltime()
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_inverse, matrix, &
0.0_dp, tmp1, flop=flop1, filter_eps=filter_eps)
0.0_dp, tmp1, flop=flop1, filter_eps=my_filter_eps)
mymat(1)%matrix => tmp1
CALL setup_arnoldi_data(my_arnoldi, mymat, max_iter=30, threshold=1.0E-3_dp, selection_crit=1, &
nval_request=2, nrestarts=2, generalized_ev=.FALSE., iram=.TRUE.)
CALL arnoldi_ev(mymat, my_arnoldi)
max_eV = REAL(get_selected_ritz_val(my_arnoldi, 2), dp)
min_eV = REAL(get_selected_ritz_val(my_arnoldi, 1), dp)
CALL deallocate_arnoldi_data(my_arnoldi)
IF (accelerator_type == 1) THEN
occ_matrix = dbcsr_get_occupation(matrix_inverse)
! 2.0 would be the correct scaling howver, we should make sure here, that we are in the convergence radius
CALL dbcsr_scale(tmp1, 1.9_dp/(min_ev+max_ev))
CALL dbcsr_scale(matrix_inverse, 1.9_dp/(min_ev+max_ev))
min_ev = min_ev*1.9_dp/(min_ev+max_ev)
! scale the matrix to get into the convergence range
CALL arnoldi_extremal(tmp1, max_eV, min_eV, threshold=my_eps_lanczos, &
max_iter=my_max_iter_lanczos, converged=arnoldi_converged)
!mymat(1)%matrix => tmp1
!CALL setup_arnoldi_data(my_arnoldi, mymat, max_iter=30, threshold=1.0E-3_dp, selection_crit=1, &
! nval_request=2, nrestarts=2, generalized_ev=.FALSE., iram=.TRUE.)
!CALL arnoldi_ev(mymat, my_arnoldi)
!max_eV = REAL(get_selected_ritz_val(my_arnoldi, 2), dp)
!min_eV = REAL(get_selected_ritz_val(my_arnoldi, 1), dp)
!CALL deallocate_arnoldi_data(my_arnoldi)
IF (unit_nr > 0) THEN
WRITE (unit_nr, *)
WRITE (unit_nr, '(T6,A,1X,L1,A,E12.3)') "Lanczos converged: ", arnoldi_converged, " threshold:", my_eps_lanczos
WRITE (unit_nr, '(T6,A,1X,E12.3,E12.3)') "Est. extremal eigenvalues:", max_eV, min_eV
WRITE (unit_nr, '(T6,A,1X,E12.3)') "Est. condition number :", max_eV/MAX(min_eV, EPSILON(min_eV))
ENDIF
! 2.0 would be the correct scaling however, we should make sure here, that we are in the convergence radius
scalingf = 1.9_dp/(max_eV+min_eV)
CALL dbcsr_scale(tmp1, scalingf)
CALL dbcsr_scale(matrix_inverse, scalingf)
min_ev = min_ev*scalingf
ENDIF
! done with the initial guess, start iterations
converged = .FALSE.
DO i = 1, 100
! tmp1 = S^-1 S
! for the convergence check
!frob_matrix_base=dbcsr_frobenius_norm(tmp1)
CALL dbcsr_add_on_diag(tmp1, -1.0_dp)
frob_matrix = dbcsr_frobenius_norm(tmp1)
CALL dbcsr_norm(tmp1, &
dbcsr_norm_maxabsnorm, norm_scalar=maxnorm_matrix)
CALL dbcsr_add_on_diag(tmp1, +1.0_dp)
! tmp2 = S^-1 S S^-1
CALL dbcsr_multiply("N", "N", 1.0_dp, tmp1, matrix_inverse, 0.0_dp, tmp2, &
flop=flop2, filter_eps=filter_eps)
flop=flop2, filter_eps=my_filter_eps)
! S^-1_{n+1} = 2 S^-1 - S^-1 S S^-1
CALL dbcsr_add(matrix_inverse, tmp2, 2.0_dp, -1.0_dp)
CALL dbcsr_filter(matrix_inverse, threshold)
CALL dbcsr_filter(matrix_inverse, my_filter_eps)
t2 = m_walltime()
occ_matrix = dbcsr_get_occupation(matrix_inverse)
! use the scalar form of the algorithm to trace the EV
min_ev = min_ev*(2.0_dp-min_ev)
IF (PRESENT(norm_convergence)) maxnorm_matrix = ABS(min_eV-1.0_dp)
IF (accelerator_type == 1) THEN
min_ev = min_ev*(2.0_dp-min_ev)
IF (PRESENT(norm_convergence)) maxnorm_matrix = ABS(min_eV-1.0_dp)
ENDIF
IF (unit_nr > 0) THEN
WRITE (unit_nr, '(T6,A,1X,I3,1X,F10.8,E12.3,F12.3,F13.3)') "Hotelling iter", i, occ_matrix, &
@ -179,18 +598,27 @@ CONTAINS
CALL m_flush(unit_nr)
ENDIF
IF (maxnorm_matrix < convergence) EXIT
IF (maxnorm_matrix < convergence) THEN
converged = .TRUE.
EXIT
ENDIF
! scale the matrix for improved convergence
min_ev = min_ev*2.0_dp/(min_ev+1.0_dp)
CALL dbcsr_scale(matrix_inverse, 2.0_dp/(min_ev+1.0_dp))
IF (accelerator_type == 1) THEN
min_ev = min_ev*2.0_dp/(min_ev+1.0_dp)
CALL dbcsr_scale(matrix_inverse, 2.0_dp/(min_ev+1.0_dp))
ENDIF
t1 = m_walltime()
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_inverse, matrix, &
0.0_dp, tmp1, flop=flop1, filter_eps=filter_eps)
0.0_dp, tmp1, flop=flop1, filter_eps=my_filter_eps)
ENDDO
IF (.NOT. converged) THEN
CPABORT("Hotelling inversion did not converge")
ENDIF
! try to symmetrize the output matrix
IF (dbcsr_get_matrix_type(matrix_inverse) == dbcsr_type_no_symmetry) THEN
CALL dbcsr_transposed(tmp2, matrix_inverse)
@ -349,8 +777,9 @@ CONTAINS
INTEGER(KIND=int_8) :: flop1, flop2, flop3, flop4, flop5
LOGICAL :: arnoldi_converged
REAL(KIND=dp) :: a, b, c, conv, d, frob_matrix, &
frob_matrix_base, max_ev, min_ev, oa, &
ob, oc, occ_matrix, od, scaling, t1, t2
frob_matrix_base, gershgorin_norm, &
max_ev, min_ev, oa, ob, oc, &
occ_matrix, od, scaling, t1, t2
TYPE(cp_logger_type), POINTER :: logger
TYPE(dbcsr_type) :: tmp1, tmp2, tmp3
@ -378,17 +807,28 @@ CONTAINS
CALL dbcsr_copy(matrix_sqrt, matrix)
! scale the matrix to get into the convergence range
CALL arnoldi_extremal(matrix_sqrt, max_ev, min_ev, threshold=eps_lanczos, &
max_iter=max_iter_lanczos, converged=arnoldi_converged)
IF (unit_nr > 0) THEN
WRITE (unit_nr, *)
WRITE (unit_nr, '(T6,A,1X,L1,A,E12.3)') "Lanczos converged: ", arnoldi_converged, " threshold:", eps_lanczos
WRITE (unit_nr, '(T6,A,1X,E12.3,E12.3)') "Est. extremal eigenvalues:", max_ev, min_ev
WRITE (unit_nr, '(T6,A,1X,E12.3)') "Est. condition number :", max_ev/MAX(min_ev, EPSILON(min_ev))
IF (order == 0) THEN
gershgorin_norm = dbcsr_gershgorin_norm(matrix_sqrt)
frob_matrix = dbcsr_frobenius_norm(matrix_sqrt)
scaling = 1.0_dp/MIN(frob_matrix, gershgorin_norm)
ELSE
! scale the matrix to get into the convergence range
CALL arnoldi_extremal(matrix_sqrt, max_ev, min_ev, threshold=eps_lanczos, &
max_iter=max_iter_lanczos, converged=arnoldi_converged)
IF (unit_nr > 0) THEN
WRITE (unit_nr, *)
WRITE (unit_nr, '(T6,A,1X,L1,A,E12.3)') "Lanczos converged: ", arnoldi_converged, " threshold:", eps_lanczos
WRITE (unit_nr, '(T6,A,1X,E12.3,E12.3)') "Est. extremal eigenvalues:", max_ev, min_ev
WRITE (unit_nr, '(T6,A,1X,E12.3)') "Est. condition number :", max_ev/MAX(min_ev, EPSILON(min_ev))
ENDIF
! conservatively assume we get a relatively large error (100*threshold_lanczos) in the estimates
! and adjust the scaling to be on the safe side
scaling = 2.0_dp/(max_ev+min_ev+100*eps_lanczos)
ENDIF
! conservatively assume we get a relatively large error (100*threshold_lanczos) in the estimates
! and adjust the scaling to be on the safe side
scaling = 2/(max_ev+min_ev+100*eps_lanczos)
CALL dbcsr_scale(matrix_sqrt, scaling)
CALL dbcsr_filter(matrix_sqrt, threshold)
@ -408,7 +848,7 @@ CONTAINS
flop4 = 0; flop5 = 0
SELECT CASE (order)
CASE (2)
CASE (0, 2)
! update the above to 0.5*(3*I-Zk*Yk)
CALL dbcsr_add_on_diag(tmp1, -2.0_dp)
CALL dbcsr_scale(tmp1, -0.5_dp)
@ -465,7 +905,7 @@ CONTAINS
! final scale
CALL dbcsr_scale(tmp1, 35.0_dp/128.0_dp)
CASE DEFAULT
CPABORT("")
CPABORT("Illegal order value")
END SELECT
! tmp2 = Yk * tmp1 = Y(k+1)

3
src/motion/integrator.F
View file

@ -215,6 +215,9 @@ CONTAINS
! atoms because of the possiblity of Langevin regions, and var_w
! for each region should depend on the temperature defined in the
! region
! RZK explains: sigma is the variance of the Wiener process associated
! with the stochastic term, sigma = m*var_w = m*(2*k_B*T*gamma*dt),
! noisy_gamma adds excessive noise that is not balanced by the damping term
ALLOCATE(var_w(nparticle))
var_w(1:nparticle) = simpar%var_w
IF (simpar%do_thermal_region) THEN

9
src/mscfg_methods.F
View file

@ -374,7 +374,7 @@ CONTAINS
INTEGER :: almo_guess_type, frz_term_type, &
method_name_id, scf_guess_type
LOGICAL :: is_crystal, is_fast_dirty
LOGICAL :: almo_scf_is_on, is_crystal, is_fast_dirty
TYPE(molecular_scf_guess_env_type), POINTER :: mscfg_env
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(section_vals_type), POINTER :: force_env_section, subsection
@ -383,6 +383,7 @@ CONTAINS
! What kind of options are we using in the loop ?
is_fast_dirty = .TRUE.
is_crystal = .FALSE.
almo_scf_is_on = .FALSE.
NULLIFY (qs_env, mscfg_env, force_env_section, subsection)
CALL force_env_get(force_env, force_env_section=force_env_section)
@ -405,6 +406,10 @@ CONTAINS
NULLIFY (subsection)
subsection => section_vals_get_subs_vals(force_env_section, "DFT%ALMO_SCF")
CALL section_vals_val_get(subsection, "ALMO_SCF_GUESS", i_val=almo_guess_type)
! check whether ALMO SCF is on
NULLIFY (subsection)
subsection => section_vals_get_subs_vals(force_env_section, "DFT%QS")
CALL section_vals_val_get(subsection, "ALMO_SCF", l_val=almo_scf_is_on)
! check SCF guess option
NULLIFY (subsection)
@ -419,7 +424,7 @@ CONTAINS
! Are we doing the loop ?
IF (scf_guess_type .EQ. molecular_guess .OR. & ! SCF guess is molecular
almo_guess_type .EQ. molecular_guess .OR. & ! ALMO SCF guess is molecular
(almo_guess_type .EQ. molecular_guess .AND. almo_scf_is_on) .OR. & ! ALMO SCF guess is molecular
frz_term_type .NE. almo_frz_none) THEN ! ALMO FRZ term is requested
do_mol_loop = .TRUE.

11
src/optimize_basis.F
View file

@ -289,14 +289,17 @@ CONTAINS
qs_env%admm_env%work_aux_aux, energy(my_id))
my_time(my_id) = m_walltime()-start_time(icalc)
IF (.NOT. para_env%ionode) THEN
f_vec = 0.0_dp; cond_vec = 0.0_dp; my_time = 0.0_dp; energy = 0.0_dp
END IF
END DO
IF (.NOT. para_env%ionode) THEN
f_vec = 0.0_dp; cond_vec = 0.0_dp; my_time = 0.0_dp; energy = 0.0_dp
END IF
DEALLOCATE (start_time)
CALL mp_sum(f_vec, para_env_top%group); CALL mp_sum(cond_vec, para_env_top%group); CALL mp_sum(my_time, para_env_top%group)
CALL mp_sum(f_vec, para_env_top%group)
CALL mp_sum(cond_vec, para_env_top%group)
CALL mp_sum(my_time, para_env_top%group)
CALL mp_sum(energy, para_env_top%group)
opt_bas%powell_param%f = 0.0_dp
DO icalc = 1, SIZE(f_vec)
icomb = MOD(icalc-1, opt_bas%ncombinations)

16
src/qs_force.F
View file

@ -310,22 +310,6 @@ CONTAINS
! Compute grid-based forces
CALL qs_ks_update_qs_env(qs_env, calculate_forces=.TRUE.)
! ALMO Code (in the spirit of the MP2 modifications below)
IF (ASSOCIATED(qs_env%almo_scf_env)) THEN
! tell qs about the energy correction
NULLIFY (energy)
CALL get_qs_env(qs_env, energy=energy)
energy%total = energy%total+energy%singles_corr
END IF
! ! ALMO Code (in the spirit of the MP2 modifications below)
! IF (ASSOCIATED(qs_env%almo_scf_env)) THEN
! ! tell qs about the energy correction
! NULLIFY (energy)
! CALL get_qs_env(qs_env, energy=energy)
! energy%total = energy%total+energy%singles_corr
! END IF
! MP2 Code
IF (ASSOCIATED(qs_env%mp2_env)) THEN
NULLIFY (matrix_p_mp2, matrix_w_mp2, rho, ks_env, energy)

47
src/qs_subsys_methods.F
View file

@ -162,11 +162,10 @@ CONTAINS
routineP = moduleN//':'//routineN
INTEGER :: arbitrary_spin, iatom, ikind, imol, &
ispin, n_ao, natom, nmol_kind, nsgf, &
nspins, z_molecule
n_ao, natom, nmol_kind, nsgf, nspins, &
z_molecule
INTEGER, DIMENSION(0:lmat, 10) :: ne_core, ne_elem, ne_explicit
INTEGER, DIMENSION(2) :: n_elec_alpha_and_beta, &
n_occ_alpha_and_beta
INTEGER, DIMENSION(2) :: n_occ_alpha_and_beta
REAL(KIND=dp) :: charge_molecule, zeff, zeff_correction
REAL(KIND=dp), DIMENSION(0:lmat, 10, 2) :: edelta
TYPE(all_potential_type), POINTER :: all_potential
@ -183,6 +182,7 @@ CONTAINS
natom = 0
! *** Initialize the molecule kind data structure ***
ARBITRARY_SPIN = 1
DO imol = 1, nmol_kind
molecule_kind => molecule_kind_set(imol)
@ -191,7 +191,6 @@ CONTAINS
!nelectron = 0
n_ao = 0
n_occ_alpha_and_beta(1:nspins) = 0
n_elec_alpha_and_beta(1:nspins) = 0
z_molecule = 0
DO iatom = 1, natom
@ -204,25 +203,30 @@ CONTAINS
gth_potential=gth_potential, &
sgp_potential=sgp_potential)
! Get the electronic state of the atom that is used
! to calculate the ATOMIC GUESS
! The spin loop is becsause ATOMIC guess can be calculated twice
! for two separate states with their own alpha-beta combinations
! This is done to break the spin symmetry of the initial wfn
! Obtain the electronic state of the atom
! The same state is used to calculate the ATOMIC GUESS
! It is great that we are consistent with ATOMIC_GUESS
CALL init_atom_electronic_state(atomic_kind=atomic_kind, &
qs_kind=qs_kind_set(ikind), &
ncalc=ne_explicit, &
ncore=ne_core, &
nelem=ne_elem, &
edelta=edelta)
DO ispin = 1, nspins
! Get the number of electrons: explicit (i.e. with orbitals) and total
! Note that it is impossible to separate alpha and beta electrons
! because the mupliplicity of the atomic state is not specified in ATOMIC GUESS
n_occ_alpha_and_beta(ispin) = n_occ_alpha_and_beta(ispin)+SUM(ne_explicit)+ &
SUM(NINT(edelta(:, :, ispin)))
n_elec_alpha_and_beta(ispin) = n_elec_alpha_and_beta(ispin)+SUM(ne_elem)
ENDDO
! If &BS section is used ATOMIC_GUESS is calculated twice
! for two separate wfns with their own alpha-beta combinations
! This is done to break the spin symmetry of the initial wfn
! For now, only alpha part of &BS is used to count electrons on
! molecules
! Get the number of explicit electrons (i.e. with orbitals)
! For now, only the total number of electrons can be obtained
! from init_atom_electronic_state
n_occ_alpha_and_beta(ARBITRARY_SPIN) = &
n_occ_alpha_and_beta(ARBITRARY_SPIN)+SUM(ne_explicit)+ &
SUM(NINT(2*edelta(:, :, ARBITRARY_SPIN)))
! We need a way to specify the number of alpha and beta electrons
! on each molecule (i.e. multiplicity is not enough)
!n_occ(ispin) = n_occ(ispin) + SUM(ne_explicit) + SUM(NINT(2*edelta(:, :, ispin)))
IF (ASSOCIATED(all_potential)) THEN
CALL get_potential(potential=all_potential, zeff=zeff, &
@ -254,13 +258,6 @@ CONTAINS
! At this point we have the number of electrons (alpha+beta) on the molecule
! as they are seen by the ATOMIC GUESS routines
! First compare two copies of the ATOMIC GUESS
! If they are different warn the user that only the ALPHA copy is currently
! assigned to the molecule
ARBITRARY_SPIN = 1
! IF( n_occ_alpha_and_beta(1).ne.n_occ_alpha_and_beta(2) ) THEN
! CPErrorMessage(cp_failure_level,routineP,"SECOND SPIN CONFIG IS IGNORED WHEN MOLECULAR STATES ARE ASSIGNED!")
! END IF
charge_molecule = REAL(z_molecule-n_occ_alpha_and_beta(ARBITRARY_SPIN), dp)
CALL set_molecule_kind(molecule_kind=molecule_kind, &
nelectron=n_occ_alpha_and_beta(ARBITRARY_SPIN), &

4
tests/QS/regtest-admm-dm/TEST_FILES
View file

@ -5,7 +5,7 @@ CH3-BP-McWeeny.inp 1 2e-13 -
CH3-BP-NONE_DM.inp 1 1e-13 -7.36784986949224
CH3-BP-NONE_DM_OT_OFF.inp 1 1.0E-14 -7.3980478715044704
CH4-BP-NONE_DM.inp 1 2e-13 -8.07630758920049
CH4-BP-NONE_DM_OT_OFF.inp 1 2e-14 -8.0771213536441095
CH4-BP-NONE_DM_OT_OFF.inp 1 2e-14 -8.07712135364430
2H2O-BLOCKED-NONE_DM.inp 1 3e-14 -34.077034488245907
H2+-BLOCKED-NONE_DM.inp 1 3e-13 -0.45795647554813002
H2+-BLOCKED-NONE_DM.inp 1 3e-13 -0.45795647554813002
#EOF

14
tests/QS/regtest-almo-1/TEST_FILES
View file

@ -1,8 +1,8 @@
almo-x.inp 11 1e-11 -137.653134569590435
almo-guess.inp 11 1e-11 -137.653034945987258
almo-scf.inp 11 1e-11 -137.652212420125409
almo-d.inp 11 9e-12 -85.963926462173887
almo-fullx.inp 11 1e-12 -137.653935641006086
almo-fullx-then-scf.inp 11
almo-then-wannier.inp 11 1e-12 -137.653115136269236
almo-x.inp 11 1e-11 -137.653122944336531
almo-guess.inp 11 1e-11 -137.653042893416284
almo-scf.inp 11 1e-11 -137.652214083887657
almo-d.inp 11 9e-12 -85.963944793676802
almo-fullx.inp 11 1e-12 -137.653899517584335
almo-fullx-then-scf.inp 11 1e-12 -137.653043387705907
almo-then-wannier.inp 11 1e-12 -137.653091046048360
#EOF

1
tests/QS/regtest-almo-1/almo-then-wannier.inp
View file

@ -17,6 +17,7 @@
EPS_FILTER 1.0E-8
ALMO_ALGORITHM DIAG
ALMO_SCF_GUESS ATOMIC
RETURN_ORTHOGONALIZED_MOS TRUE
&ALMO_OPTIMIZER_DIIS
MAX_ITER 10

5
tests/QS/regtest-almo-1/almo-x.inp
View file

@ -17,6 +17,11 @@
EPS_FILTER 1.0E-8
ALMO_ALGORITHM DIAG
ALMO_SCF_GUESS ATOMIC
MO_OVERLAP_INV_ALG LS_TAYLOR
&MATRIX_ITERATE
EPS_TARGET_FACTOR 1
&END
&ALMO_OPTIMIZER_DIIS
MAX_ITER 10

3
tests/QS/regtest-almo-2/FH-chain.inp
View file

@ -17,6 +17,7 @@
EPS_FILTER 1.0E-8
ALMO_ALGORITHM DIAG
ALMO_SCF_GUESS ATOMIC
RETURN_ORTHOGONALIZED_MOS F
&ALMO_OPTIMIZER_DIIS
MAX_ITER 30
@ -33,7 +34,7 @@
CONJUGATOR FLETCHER
LIN_SEARCH_EPS_ERROR 0.05
LIN_SEARCH_STEP_SIZE_GUESS 0.1
MAX_ITER_OUTER_LOOP 0
MAX_ITER_OUTER_LOOP 2
&END XALMO_OPTIMIZER_PCG
&END ALMO_SCF

10
tests/QS/regtest-almo-2/TEST_FILES
View file

@ -1,6 +1,8 @@
almo-fullx.inp 11 1e-11 -137.653972983204994
almo-no-deloc.inp 11 9e-12 -137.557304851561014
almo-fullx.inp 11 1e-11 -137.653869649992146
almo-no-deloc.inp 11 9e-12 -137.557304848837759
FH-chain.inp 11 3e-10 -98.108319997143496
ion-pair.inp 11 1e-13 -115.006279154085647
LiF.inp 11 3e-12 -127.53931597060651
ion-pair.inp 11 1e-13 -115.006279649267199
LiF.inp 11 3e-12 -127.539319485641016
zero-electron-frag.inp 11 1e-08 -98.728605260047331
matrix-iterate.inp 11 1e-09 -127.379124535900814
#EOF

117
tests/QS/regtest-almo-2/matrix-iterate.inp Normal file
View file

@ -0,0 +1,117 @@
&GLOBAL
PROJECT iterate
RUN_TYPE ENERGY
PRINT_LEVEL LOW
&END GLOBAL
&FORCE_EVAL
METHOD QS
&DFT
BASIS_SET_FILE_NAME GTH_BASIS_SETS
POTENTIAL_FILE_NAME POTENTIAL
&MGRID
CUTOFF 200
NGRIDS 4
&END MGRID
&QS
ALMO_SCF T
EPS_DEFAULT 1.0E-9
&END QS
&ALMO_SCF
EPS_FILTER 1.0E-8
ALMO_ALGORITHM DIAG
DELOCALIZE_METHOD XALMO_SCF
XALMO_R_CUTOFF_FACTOR 0.6
ALMO_SCF_GUESS ATOMIC
&MATRIX_ITERATE
EPS_TARGET_FACTOR 100.
EPS_LANCZOS 1.0E-3
MAX_ITER_LANCZOS 64
ORDER_LANCZOS 3
&END MATRIX_ITERATE
&ALMO_OPTIMIZER_DIIS
MAX_ITER 100
N_DIIS 5
EPS_ERROR 1.0E-5
&END ALMO_OPTIMIZER_DIIS
&XALMO_OPTIMIZER_PCG
MAX_ITER 30
EPS_ERROR 1.0E-7
CONJUGATOR HESTENES_STIEFEL
PRECONDITIONER DEFAULT
LIN_SEARCH_EPS_ERROR 0.05
LIN_SEARCH_STEP_SIZE_GUESS 0.1
MAX_ITER_OUTER_LOOP 2
&END XALMO_OPTIMIZER_PCG
&END ALMO_SCF
&XC
&XC_FUNCTIONAL PBE
&END XC_FUNCTIONAL
&END XC
&END DFT
&SUBSYS
&CELL
ABC 4.0351 4.0351 4.0351
MULTIPLE_UNIT_CELL 1 1 1
&END CELL
&TOPOLOGY
MULTIPLE_UNIT_CELL 1 1 1
&END
&COORD
SCALED
Li 0.0 0.0 0.0 Li-plus
Li 0.5 0.5 0.0 Li-plus
Li 0.5 0.0 0.5 Li-plus
Li 0.0 0.5 0.5 Li-plus
F 0.0 0.0 0.5 F-minus
F 0.0 0.5 0.0 F-minus
F 0.5 0.0 0.0 F-minus
F 0.5 0.5 0.5 F-minus
&END COORD
&KIND Li
BASIS_SET DZVP-GTH
POTENTIAL GTH-BLYP-q3
&BS
&ALPHA
NEL -1
L 0
N 2
&END
! BETA FUNCTION SHOULD BE THE SAME
! TO AVOID WARNINGS
&BETA
NEL -1
L 0
N 2
&END
&END
&END KIND
&KIND F
BASIS_SET DZVP-GTH
POTENTIAL GTH-BLYP-q7
&BS
&ALPHA
NEL +1
L 1
N 2
&END
! BETA FUNCTION SHOULD BE THE SAME
! TO AVOID WARNINGS
&BETA
NEL +1
L 1
N 2
&END
&END
&END KIND
&END SUBSYS
&END FORCE_EVAL

112
tests/QS/regtest-almo-2/zero-electron-frag.inp Normal file
View file

@ -0,0 +1,112 @@
&GLOBAL
PROJECT LiF-chain
RUN_TYPE ENERGY
PRINT_LEVEL LOW
&END GLOBAL
&FORCE_EVAL
METHOD QS
&DFT
POTENTIAL_FILE_NAME GTH_POTENTIALS
BASIS_SET_FILE_NAME GTH_BASIS_SETS
&QS
ALMO_SCF T
EPS_DEFAULT 1.0E-8 ! 1.0E-12
&END QS
&ALMO_SCF
EPS_FILTER 1.0E-8
ALMO_ALGORITHM DIAG
ALMO_SCF_GUESS ATOMIC
&ALMO_OPTIMIZER_DIIS
MAX_ITER 30
EPS_ERROR 5.0E-4
N_DIIS 7
&END ALMO_OPTIMIZER_DIIS
DELOCALIZE_METHOD XALMO_SCF
XALMO_R_CUTOFF_FACTOR 1.2
&XALMO_OPTIMIZER_PCG
MAX_ITER 100
EPS_ERROR 5.0E-4
CONJUGATOR FLETCHER
LIN_SEARCH_EPS_ERROR 0.05
LIN_SEARCH_STEP_SIZE_GUESS 0.1
MAX_ITER_OUTER_LOOP 0
&END XALMO_OPTIMIZER_PCG
&END ALMO_SCF
&MGRID
CUTOFF 200 ! 320
NGRIDS 5
&END MGRID
&XC
&XC_FUNCTIONAL BLYP
&END XC_FUNCTIONAL
&END XC
&END DFT
&SUBSYS
&CELL
ABC 5.0000000000 5.0000000000 10.0000000000
MULTIPLE_UNIT_CELL 1 1 1
&END CELL
&TOPOLOGY
MULTIPLE_UNIT_CELL 1 1 1
&END
&COORD
! atomic decomposition
H 0.0000000000 0.0000000000 0.0000000000 H1
F 0.0000000000 0.0000000000 1.0000000000 F1
H 0.0000000000 0.0000000000 3.0000000000 H2
F 0.0000000000 1.0000000000 3.0000000000 F2
H 0.0000000000 0.0000000000 6.0000000000 H3
F 0.0000000000 0.0000000000 5.0000000000 F3
H 0.0000000000 1.0000000000 8.0000000000 H4
F 0.0000000000 0.0000000000 8.0000000000 F4
! molecular decomposition
!H 0.0000000000 0.0000000000 0.0000000000 HF1
!F 0.0000000000 0.0000000000 1.0000000000 HF1
!H 0.0000000000 0.0000000000 3.0000000000 HF2
!F 0.0000000000 1.0000000000 3.0000000000 HF2
!H 0.0000000000 0.0000000000 6.0000000000 HF3
!F 0.0000000000 0.0000000000 5.0000000000 HF3
!H 0.0000000000 1.0000000000 8.0000000000 HF4
!F 0.0000000000 0.0000000000 8.0000000000 HF4
&END COORD
&KIND H
BASIS_SET DZVP-GTH
POTENTIAL GTH-BLYP-q1
&BS
&ALPHA
NEL -1
L 0
N 1
&END
&BETA
NEL -1
L 0
N 1
&END
&END
&END KIND
&KIND F
BASIS_SET DZVP-GTH
POTENTIAL GTH-BLYP-q7
&BS
&ALPHA
NEL +1
L 1
N 2
&END
&BETA
NEL +1
L 1
N 2
&END
&END
&END KIND
&END SUBSYS
&END FORCE_EVAL

2
tests/QS/regtest-almo-eda/TEST_FILES
View file

@ -1,2 +1,2 @@
almo-eda-x.inp 11 1e-08 -85.963728515581991
almo-eda-x.inp 11 1e-08 -85.963743073641155
#EOF

12
tests/QS/regtest-almo-md/TEST_FILES
View file

@ -1,7 +1,7 @@
almo-md-full-scf.inp 11
almo-md-full-x-then-scf.inp 11
almo-md-no-aspc.inp 11
xalmo-scf-md.inp 11
almo-md.inp 11
almo-md-wannier.inp 11
almo-md-full-scf.inp 11 1.e-12 -136.705845206653748
almo-md-full-x-then-scf.inp 11 1.e-12 -136.986095053499298
almo-md-no-aspc.inp 11 1.e-12 -136.473484271103160
xalmo-scf-md.inp 11 1.e-12 -136.705599282831827
almo-md.inp 11 1.e-12 -136.622582483465919
almo-md-wannier.inp 11 1.e-12 -136.473478972490227
#EOF

4
tests/QS/regtest-almo-md/almo-md-wannier.inp
View file

@ -36,7 +36,9 @@
DELOCALIZE_METHOD NONE
XALMO_R_CUTOFF_FACTOR 1.4
RETURN_ORTHOGONALIZED_MOS TRUE
&XALMO_OPTIMIZER_PCG
MAX_ITER 100
EPS_ERROR 5.0E-4

147
tests/QS/regtest-almo-strong/BNH2-sheet.inp Normal file
View file

@ -0,0 +1,147 @@
&GLOBAL
PROJECT BNH
RUN_TYPE ENERGY
PRINT_LEVEL LOW
!TRACE TRUE
&END GLOBAL
&FORCE_EVAL
METHOD QS
STRESS_TENSOR ANALYTICAL
&DFT
POTENTIAL_FILE_NAME GTH_POTENTIALS
BASIS_SET_FILE_NAME GTH_BASIS_SETS
&QS
ALMO_SCF T
EPS_DEFAULT 1.E-8 !1.0E-14
&END QS
&ALMO_SCF
EPS_FILTER 1.0E-8 !1.0E-12
ALMO_ALGORITHM DIAG
ALMO_SCF_GUESS ATOMIC
&ALMO_OPTIMIZER_DIIS
MAX_ITER 100
EPS_ERROR 1.0E-3 !1.0E-6
N_DIIS 5
&END ALMO_OPTIMIZER_DIIS
DELOCALIZE_METHOD XALMO_SCF
XALMO_R_CUTOFF_FACTOR 0.6
&XALMO_OPTIMIZER_PCG
MAX_ITER 100
EPS_ERROR 1.0E-3 !1.0E-6
CONJUGATOR FLETCHER_REEVES
LIN_SEARCH_EPS_ERROR 0.01
LIN_SEARCH_STEP_SIZE_GUESS 0.2
MAX_ITER_OUTER_LOOP 0
&END XALMO_OPTIMIZER_PCG
&END ALMO_SCF
&MGRID
CUTOFF 200 !600
NGRIDS 5
&END MGRID
&XC
&XC_FUNCTIONAL PBE
&END XC_FUNCTIONAL
&END XC
&END DFT
&SUBSYS
&CELL
ABC 2.647 2.647 5.0
ALPHA_BETA_GAMMA 90.0 90.0 120.0
MULTIPLE_UNIT_CELL 1 1 1
&END CELL
&TOPOLOGY
!&GENERATE
! BONDLENGTH_MAX 5.2
!&END GENERATE
MULTIPLE_UNIT_CELL 1 1 1
&END
&COORD
B 0.0000 0.0000 0.1936 B1
Hm 0.0000 0.0000 1.3944 B1
N 0.0000 1.5282 -0.2704 N1
Hp 0.0000 1.5282 -1.3173 N1
&END COORD
&KIND Hp
ELEMENT H
BASIS_SET SZV-GTH
POTENTIAL GTH-BLYP-q1
!&BS
! &ALPHA
! NEL -1
! L 0
! N 1
! &END
! ! BETA FUNCTION SHOULD BE THE SAME
! ! TO AVOID WARNINGS
! &BETA
! NEL -1
! L 0
! N 1
! &END
!&END
&END KIND
&KIND Hm
ELEMENT H
BASIS_SET SZV-GTH
POTENTIAL GTH-BLYP-q1
!&BS
! &ALPHA
! NEL 1
! L 0
! N 1
! &END
! ! BETA FUNCTION SHOULD BE THE SAME
! ! TO AVOID WARNINGS
! &BETA
! NEL 1
! L 0
! N 1
! &END
!&END
&END KIND
&KIND N
BASIS_SET SZV-GTH
POTENTIAL GTH-BLYP-q5
&BS
&ALPHA
NEL 2
L 1
N 2
&END
! BETA FUNCTION SHOULD BE THE SAME
! TO AVOID WARNINGS
&BETA
NEL 2
L 1
N 2
&END
&END
&END KIND
&KIND B
BASIS_SET SZV-GTH
POTENTIAL GTH-BLYP-q3
&BS
&ALPHA
NEL -1 -1
L 0 1
N 2 2
&END
! BETA FUNCTION SHOULD BE THE SAME
! TO AVOID WARNINGS
&BETA
NEL -1 -1
L 0 1
N 2 2
&END
&END
&END KIND
&END SUBSYS
&END FORCE_EVAL

113
tests/QS/regtest-almo-strong/FH-chain.strong.inp Normal file
View file

@ -0,0 +1,113 @@
&GLOBAL
PROJECT LiF-chain
RUN_TYPE ENERGY
PRINT_LEVEL LOW
&END GLOBAL
&FORCE_EVAL
METHOD QS
&DFT
POTENTIAL_FILE_NAME GTH_POTENTIALS
BASIS_SET_FILE_NAME GTH_BASIS_SETS
&QS
ALMO_SCF T
EPS_DEFAULT 1.0E-8 ! 1.0E-12
&END QS
&ALMO_SCF
EPS_FILTER 1.0E-8
ALMO_ALGORITHM SKIP
ALMO_SCF_GUESS ATOMIC
XALMO_TRIAL_WF SIMPLE
&ALMO_OPTIMIZER_DIIS
MAX_ITER 30
EPS_ERROR 5.0E-4
N_DIIS 7
&END ALMO_OPTIMIZER_DIIS
DELOCALIZE_METHOD XALMO_SCF
XALMO_R_CUTOFF_FACTOR 1.2
&XALMO_OPTIMIZER_PCG
MAX_ITER 50
EPS_ERROR 5.0E-4
CONJUGATOR FLETCHER
LIN_SEARCH_EPS_ERROR 0.05
LIN_SEARCH_STEP_SIZE_GUESS 0.1
MAX_ITER_OUTER_LOOP 2
&END XALMO_OPTIMIZER_PCG
&END ALMO_SCF
&MGRID
CUTOFF 200 ! 320
NGRIDS 5
&END MGRID
&XC
&XC_FUNCTIONAL BLYP
&END XC_FUNCTIONAL
&END XC
&END DFT
&SUBSYS
&CELL
ABC 3.0000000000 3.0000000000 10.0000000000
MULTIPLE_UNIT_CELL 1 1 1
&END CELL
&TOPOLOGY
MULTIPLE_UNIT_CELL 1 1 1
&END
&COORD
! atomic decomposition
H 0.0000000000 0.0000000000 0.0000000000 H1
F 0.0000000000 0.0000000000 1.0000000000 F1
H 0.0000000000 0.0000000000 3.0000000000 H2
F 0.0000000000 1.0000000000 3.0000000000 F2
H 0.0000000000 0.0000000000 6.0000000000 H3
F 0.0000000000 0.0000000000 5.0000000000 F3
H 0.0000000000 1.0000000000 8.0000000000 H4
F 0.0000000000 0.0000000000 8.0000000000 F4
! molecular decomposition
!H 0.0000000000 0.0000000000 0.0000000000 HF1
!F 0.0000000000 0.0000000000 1.0000000000 HF1
!H 0.0000000000 0.0000000000 3.0000000000 HF2
!F 0.0000000000 1.0000000000 3.0000000000 HF2
!H 0.0000000000 0.0000000000 6.0000000000 HF3
!F 0.0000000000 0.0000000000 5.0000000000 HF3
!H 0.0000000000 1.0000000000 8.0000000000 HF4
!F 0.0000000000 0.0000000000 8.0000000000 HF4
&END COORD
&KIND H
BASIS_SET DZVP-GTH
POTENTIAL GTH-BLYP-q1
&BS
&ALPHA
NEL -1
L 0
N 1
&END
&BETA
NEL -1
L 0
N 1
&END
&END
&END KIND
&KIND F
BASIS_SET DZVP-GTH
POTENTIAL GTH-BLYP-q7
&BS
&ALPHA
NEL +1
L 1
N 2
&END
&BETA
NEL +1
L 1
N 2
&END
&END
&END KIND
&END SUBSYS
&END FORCE_EVAL

5
tests/QS/regtest-almo-strong/TEST_FILES Normal file
View file

@ -0,0 +1,5 @@
atomic-water.inp 11 1e-08 -33.939454719891792
bn.inp 11 1e-10 -50.599959228858097
BNH2-sheet.inp 11 1e-08 -13.641848381243719
FH-chain.strong.inp 11 1e-08 -98.638850796080803
#EOF

132
tests/QS/regtest-almo-strong/atomic-water.inp Normal file
View file

@ -0,0 +1,132 @@
&GLOBAL
PROJECT md-atomic-part
RUN_TYPE MD
PRINT_LEVEL LOW
&END GLOBAL
&MOTION
&MD
ENSEMBLE NVE
STEPS 4
TIMESTEP 0.5
TEMPERATURE 298
&END MD
&END MOTION
&FORCE_EVAL
METHOD QS
&DFT
BASIS_SET_FILE_NAME GTH_BASIS_SETS
POTENTIAL_FILE_NAME POTENTIAL
&MGRID
CUTOFF 400
&END MGRID
&QS
ALMO_SCF T
EPS_DEFAULT 1.0E-8 ! 1.0E-12
WF_INTERPOLATION PS
&END QS
&ALMO_SCF
EPS_FILTER 1.0E-8 ! 1.0E-11
ALMO_ALGORITHM SKIP
ALMO_SCF_GUESS ATOMIC
MO_OVERLAP_INV_ALG DENSE_CHOLESKY
XALMO_TRIAL_WF SIMPLE
DELOCALIZE_METHOD XALMO_SCF
XALMO_R_CUTOFF_FACTOR 1.3
XALMO_EXTRAPOLATION_ORDER 4
&MATRIX_ITERATE
EPS_TARGET_FACTOR 10.
!EPS_LANCZOS 1.0E-3
!MAX_ITER_LANCZOS 64
!ORDER_LANCZOS 3
&END MATRIX_ITERATE
&XALMO_OPTIMIZER_PCG
MAX_ITER 15
EPS_ERROR 5.0E-5
CONJUGATOR FLETCHER_REEVES
LIN_SEARCH_EPS_ERROR 0.05
LIN_SEARCH_STEP_SIZE_GUESS 0.1
MAX_ITER_OUTER_LOOP 50
PRECOND_FILTER_THRESHOLD 0.008
&END XALMO_OPTIMIZER_PCG
&END ALMO_SCF
&XC
&XC_FUNCTIONAL BLYP
&END XC_FUNCTIONAL
&END XC
&END DFT
&SUBSYS
&KIND H
BASIS_SET SZV-GTH ! TZV2P-GTH
POTENTIAL GTH-BLYP-q1
&BS
&ALPHA
NEL -1
L 0
N 1
&END
! BETA FUNCTION SHOULD BE THE SAME TO AVOID WARNINGS
&BETA
NEL -1
L 0
N 1
&END
&END
&END KIND
&KIND O
BASIS_SET SZV-GTH ! TZV2P-GTH
POTENTIAL GTH-BLYP-q6
&BS
&ALPHA
NEL +2
L 1
N 2
&END
! BETA FUNCTION SHOULD BE THE SAME TO AVOID WARNINGS
&BETA
NEL +2
L 1
N 2
&END
&END
&END KIND
&CELL
ABC 3.00000000 3.00000000 5.00000000
MULTIPLE_UNIT_CELL 1 1 1
&END CELL
&TOPOLOGY
&GENERATE # fragments contain single atoms
BONDLENGTH_MAX 1.0
BONDPARM COVALENT
BONDPARM_FACTOR 0.3
&END GENERATE
MULTIPLE_UNIT_CELL 1 1 1
&END
&COORD
O 1.5285821781 1.7062403635 3.9140529453
H 1.6296675630 1.2721242389 4.7873221686
H 2.0421991683 2.5306450214 4.0774366182
O 1.6497475797 0.9866814243 1.1082336734
H 1.5853337875 1.1818653541 2.0640644445
H 0.6898667621 0.8787694446 0.8082114051
&END COORD
&END SUBSYS
&END FORCE_EVAL

114
tests/QS/regtest-almo-strong/bn.inp Normal file
View file

@ -0,0 +1,114 @@
&GLOBAL
PROJECT BN_neg
RUN_TYPE ENERGY
PRINT_LEVEL LOW
&END GLOBAL
&FORCE_EVAL
METHOD QS
&DFT
BASIS_SET_FILE_NAME GTH_BASIS_SETS
POTENTIAL_FILE_NAME GTH_POTENTIALS
&QS
METHOD GPW
ALMO_SCF T
EPS_DEFAULT 1.0E-10 ! 1.0E-11
&END QS
&MGRID
CUTOFF 200 ! 600
NGRIDS 4
&END MGRID
&ALMO_SCF
EPS_FILTER 1.0E-9 ! 1.0E-10
ALMO_ALGORITHM SKIP
ALMO_SCF_GUESS ATOMIC
MO_OVERLAP_INV_ALG DENSE_CHOLESKY
DELOCALIZE_METHOD XALMO_SCF
XALMO_R_CUTOFF_FACTOR 0.6
XALMO_TRIAL_WF SIMPLE
RETURN_ORTHOGONALIZED_MOS FALSE
&XALMO_OPTIMIZER_PCG
MAX_ITER 15
EPS_ERROR 1.0E-5
CONJUGATOR DAI_YUAN
LIN_SEARCH_EPS_ERROR 0.01
LIN_SEARCH_STEP_SIZE_GUESS 0.05
MAX_ITER_OUTER_LOOP 10
PRECOND_FILTER_THRESHOLD 0.05
&END XALMO_OPTIMIZER_PCG
&END ALMO_SCF
&XC
&XC_FUNCTIONAL BLYP
&END XC_FUNCTIONAL
&END XC
&END DFT
&SUBSYS
&CELL
ABC 3.66 3.66 3.66
MULTIPLE_UNIT_CELL 1 1 1
&END CELL
&COORD
SCALED T
B 0.0000000000 0.0000000000 0.0000000000
B 0.50000000 0.50000000 0.0000000000
B 0.50000000 0.0000000000 0.50000000
B 0.0000000000 0.50000000 0.50000000
N 0.25000000 0.25000000 0.25000000
N 0.75000000 0.75000000 0.25000000
N 0.75000000 0.25000000 0.75000000
N 0.25000000 0.75000000 0.75000000
&END COORD
&KIND B
BASIS_SET SZV-GTH-q3 ! DZVP-GTH-q3
POTENTIAL GTH-BLYP-q3
&BS
&ALPHA
NEL -1 -2
L 1 0
N 2 2
&END
! BETA FUNCTION SHOULD BE THE SAME
! TO AVOID WARNINGS
&BETA
NEL -1 -2
L 1 0
N 2 2
&END
&END
&END KIND
&KIND N
BASIS_SET SZV-GTH-q5 ! DZVP-GTH-q5
POTENTIAL GTH-BLYP-q5
&BS
&ALPHA
NEL +3
L 1
N 2
&END
! BETA FUNCTION SHOULD BE THE SAME
! TO AVOID WARNINGS
&BETA
NEL +3
L 1
N 2
&END
&END
&END KIND
&END SUBSYS
&END FORCE_EVAL

2
tests/QS/regtest-dm-ls-scf-1/TEST_FILES
View file

@ -1,7 +1,7 @@
H2-big-1.inp 11 2e-13 -27.808422055041142
H2-big-2.inp 11 1e-13 -27.808422055528268
H2-big-3.inp 11 7e-13 -27.808586398678841
H2-big-4.inp 11 6e-13 -27.809113704800708
H2-big-4.inp 11 6e-13 -27.808505168773813
H2-big-5.inp 11 2e-13 -18.233645147453636
H2-big-6.inp 11 2e-13 -18.233645147453636
H2-big-7.inp 11 2e-13 -27.808422055041142

2
tests/QS/regtest-kg/TEST_FILES
View file

@ -22,7 +22,7 @@ H2_H2O_lsks.inp 11 1e-10 -
H2_H2O_ec.inp 66 1e-10 -18.4641086282
H2_H2O_ecprim.inp 66 1e-10 -18.4071025433
2H2O_ecmao.inp 66 1e-10 -34.0838901711
2H2O_ecmao2.inp 66 1e-08 -34.5024811487
2H2O_ecmao2.inp 66 1e-08 -34.5024586832
H2-none.inp 11 7e-07 -3.359680469888914
H2_H2O-vdW.inp 11 1e-10 -18.149003390914348
H2_H2O-lri.inp 72 3e-02 0.00011458

4
tests/QS/regtest-rma/TEST_FILES
View file

@ -1,7 +1,7 @@
H2O-32-dftb-ls-2_mult.inp 11 1e-12 -32.574187310759356
H2O-32-dftb-ls-2.inp 11 1e-12 -32.574187310759356
H2O-OT-ASPC-1.inp 1 4e-14 -17.13993294772316
H2O-OT-ASPC-1_clusters.inp 1 4e-14 -17.13993294772316
H2O-OT-ASPC-1.inp 1 4e-14 -17.13993294716181
H2O-OT-ASPC-1_clusters.inp 1 4e-14 -17.13993294716182
H2-big-nimages.inp 11 2e-13 -27.808422055041138
H2-big-nimages_clusters.inp 11 2e-13 -27.808422055041142
H2O_grad_gpw.inp 11 7e-11 -17.082584774463687

1
tests/TEST_DIRS
View file

@ -57,6 +57,7 @@ QS/regtest-admm-4 libint
ATOM/regtest-1
QS/regtest-gw-ic-model
QS/regtest-almo-md
QS/regtest-almo-strong
QS/regtest-pao-1
QS/regtest-ri-mp2 libint
QS/regtest-almo-2