mirror of
https://github.com/cp2k/cp2k.git
synced 2026-07-29 06:35:28 -04:00
PAO: Add equivariant parametrization
This commit is contained in:
Ole Schütt
Ole Schütt
parent
1a8a8ac7cf
commit
a9db994763
20 changed files with 1354 additions and 487 deletions
|
|
@ -46,6 +46,7 @@ MODULE pao_input
|
|||
pao_fock_param = 102, &
|
||||
pao_exp_param = 103, &
|
||||
pao_gth_param = 104, &
|
||||
pao_equi_param = 105, &
|
||||
pao_opt_cg = 301, &
|
||||
pao_opt_bfgs = 302, &
|
||||
pao_ml_gp = 401, &
|
||||
|
|
@ -89,83 +90,83 @@ CONTAINS
|
|||
! parse input and print
|
||||
|
||||
CALL section_vals_val_get(pao_section, "EPS_PAO", r_val=pao%eps_pao)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "EPS_PAO", pao%eps_pao
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "EPS_PAO", pao%eps_pao
|
||||
|
||||
CALL section_vals_val_get(pao_section, "MIXING", r_val=pao%mixing)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "MIXING", pao%mixing
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "MIXING", pao%mixing
|
||||
|
||||
CALL section_vals_val_get(pao_section, "MAX_PAO", i_val=pao%max_pao)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,I10)") " PAO|", "MAX_PAO", pao%max_pao
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,I11)") " PAO|", "MAX_PAO", pao%max_pao
|
||||
|
||||
CALL section_vals_val_get(pao_section, "MAX_CYCLES", i_val=pao%max_cycles)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,I10)") " PAO|", "MAX_CYCLES", pao%max_cycles
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,I11)") " PAO|", "MAX_CYCLES", pao%max_cycles
|
||||
|
||||
CALL section_vals_val_get(pao_section, "PARAMETERIZATION", i_val=pao%parameterization)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,A10)") " PAO|", "PARAMETERIZATION", id2str(pao%parameterization)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,A11)") " PAO|", "PARAMETERIZATION", id2str(pao%parameterization)
|
||||
|
||||
CALL section_vals_val_get(pao_section, "PRECONDITION", l_val=pao%precondition)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,L10)") " PAO|", "PRECONDITION", pao%precondition
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,L11)") " PAO|", "PRECONDITION", pao%precondition
|
||||
|
||||
CALL section_vals_val_get(pao_section, "REGULARIZATION", r_val=pao%regularization)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "REGULARIZATION", pao%regularization
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "REGULARIZATION", pao%regularization
|
||||
IF (pao%regularization < 0.0_dp) CPABORT("PAO: REGULARIZATION < 0")
|
||||
|
||||
CALL section_vals_val_get(input, "DFT%QS%EPS_DEFAULT", r_val=pao%eps_pgf) ! default value
|
||||
CALL section_vals_val_get(pao_section, "EPS_PGF", n_rep_val=n_rep)
|
||||
IF (n_rep /= 0) CALL section_vals_val_get(pao_section, "EPS_PGF", r_val=pao%eps_pgf)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "EPS_PGF", pao%eps_pgf
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "EPS_PGF", pao%eps_pgf
|
||||
IF (pao%eps_pgf < 0.0_dp) CPABORT("PAO: EPS_PGF < 0")
|
||||
|
||||
CALL section_vals_val_get(pao_section, "PENALTY_DISTANCE", r_val=pao%penalty_dist)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "PENALTY_DISTANCE", pao%penalty_dist
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "PENALTY_DISTANCE", pao%penalty_dist
|
||||
IF (pao%penalty_dist < 0.0_dp) CPABORT("PAO: PENALTY_DISTANCE < 0")
|
||||
|
||||
CALL section_vals_val_get(pao_section, "PENALTY_STRENGTH", r_val=pao%penalty_strength)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "PENALTY_STRENGTH", pao%penalty_strength
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "PENALTY_STRENGTH", pao%penalty_strength
|
||||
IF (pao%penalty_strength < 0.0_dp) CPABORT("PAO: PENALTY_STRENGTH < 0")
|
||||
|
||||
CALL section_vals_val_get(pao_section, "LINPOT_PRECONDITION_DELTA", r_val=pao%linpot_precon_delta)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "LINPOT_PRECONDITION_DELTA", pao%linpot_precon_delta
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "LINPOT_PRECONDITION_DELTA", pao%linpot_precon_delta
|
||||
IF (pao%linpot_precon_delta < 0.0_dp) CPABORT("PAO: LINPOT_PRECONDITION_DELTA < 0")
|
||||
|
||||
CALL section_vals_val_get(pao_section, "LINPOT_INITGUESS_DELTA", r_val=pao%linpot_init_delta)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "LINPOT_INITGUESS_DELT", pao%linpot_init_delta
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "LINPOT_INITGUESS_DELT", pao%linpot_init_delta
|
||||
IF (pao%linpot_init_delta < 0.0_dp) CPABORT("PAO: LINPOT_INITGUESS_DELTA < 0")
|
||||
|
||||
CALL section_vals_val_get(pao_section, "LINPOT_REGULARIZATION_DELTA", r_val=pao%linpot_regu_delta)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "LINPOT_REGULARIZATION_DELTA", pao%linpot_regu_delta
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "LINPOT_REGULARIZATION_DELTA", pao%linpot_regu_delta
|
||||
IF (pao%linpot_regu_delta < 0.0_dp) CPABORT("PAO: LINPOT_REGULARIZATION_DELTA < 0")
|
||||
|
||||
CALL section_vals_val_get(pao_section, "LINPOT_REGULARIZATION_STRENGTH", r_val=pao%linpot_regu_strength)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "LINPOT_REGULARIZATION_STRENGTH", pao%linpot_regu_strength
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "LINPOT_REGULARIZATION_STRENGTH", pao%linpot_regu_strength
|
||||
IF (pao%linpot_regu_strength < 0.0_dp) CPABORT("PAO: LINPOT_REGULARIZATION_STRENGTH < 0")
|
||||
|
||||
CALL section_vals_val_get(pao_section, "OPTIMIZER", i_val=pao%optimizer)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,A10)") " PAO|", "OPTIMIZER", id2str(pao%optimizer)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,A11)") " PAO|", "OPTIMIZER", id2str(pao%optimizer)
|
||||
|
||||
CALL section_vals_val_get(pao_section, "CG_INIT_STEPS", i_val=pao%cg_init_steps)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,I10)") " PAO|", "CG_INIT_STEPS", pao%cg_init_steps
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,I11)") " PAO|", "CG_INIT_STEPS", pao%cg_init_steps
|
||||
IF (pao%cg_init_steps < 1) CPABORT("PAO: CG_INIT_STEPS < 1")
|
||||
|
||||
CALL section_vals_val_get(pao_section, "CG_RESET_LIMIT", r_val=pao%cg_reset_limit)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "CG_RESET_LIMIT", pao%cg_reset_limit
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "CG_RESET_LIMIT", pao%cg_reset_limit
|
||||
IF (pao%cg_reset_limit < 0.0_dp) CPABORT("PAO: CG_RESET_LIMIT < 0")
|
||||
|
||||
CALL section_vals_val_get(pao_section, "CHECK_UNITARY_TOL", r_val=pao%check_unitary_tol)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "CHECK_UNITARY_TOL", pao%check_unitary_tol
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "CHECK_UNITARY_TOL", pao%check_unitary_tol
|
||||
|
||||
CALL section_vals_val_get(pao_section, "CHECK_GRADIENT_TOL", r_val=pao%check_grad_tol)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "CHECK_GRADIENT_TOL", pao%check_grad_tol
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "CHECK_GRADIENT_TOL", pao%check_grad_tol
|
||||
|
||||
CALL section_vals_val_get(pao_section, "NUM_GRADIENT_ORDER", i_val=pao%num_grad_order)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,I10)") " PAO|", "NUM_GRADIENT_ORDER", pao%num_grad_order
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,I11)") " PAO|", "NUM_GRADIENT_ORDER", pao%num_grad_order
|
||||
|
||||
CALL section_vals_val_get(pao_section, "NUM_GRADIENT_EPS", r_val=pao%num_grad_eps)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "NUM_GRADIENT_EPS", pao%num_grad_eps
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "NUM_GRADIENT_EPS", pao%num_grad_eps
|
||||
IF (pao%num_grad_eps < 0.0_dp) CPABORT("PAO: NUM_GRADIENT_EPS < 0")
|
||||
|
||||
CALL section_vals_val_get(pao_section, "PRINT%RESTART%WRITE_CYCLES", i_val=pao%write_cycles)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,I10)") " PAO|", "PRINT%RESTART%WRITE_CYCLES", pao%write_cycles
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,I11)") " PAO|", "PRINT%RESTART%WRITE_CYCLES", pao%write_cycles
|
||||
|
||||
CALL section_vals_val_get(pao_section, "RESTART_FILE", c_val=pao%restart_file)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,A)") " PAO|", "RESTART_FILE ", TRIM(pao%restart_file)
|
||||
|
|
@ -174,22 +175,22 @@ CONTAINS
|
|||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,A)") " PAO|", "PREOPT_DM_FILE ", TRIM(pao%preopt_dm_file)
|
||||
|
||||
CALL section_vals_val_get(pao_section, "MACHINE_LEARNING%METHOD", i_val=pao%ml_method)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,A10)") " PAO|", "MACHINE_LEARNING%METHOD", id2str(pao%ml_method)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,A11)") " PAO|", "MACHINE_LEARNING%METHOD", id2str(pao%ml_method)
|
||||
|
||||
CALL section_vals_val_get(pao_section, "MACHINE_LEARNING%PRIOR", i_val=pao%ml_prior)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,A10)") " PAO|", "MACHINE_LEARNING%PRIOR", id2str(pao%ml_prior)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,A11)") " PAO|", "MACHINE_LEARNING%PRIOR", id2str(pao%ml_prior)
|
||||
|
||||
CALL section_vals_val_get(pao_section, "MACHINE_LEARNING%DESCRIPTOR", i_val=pao%ml_descriptor)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,A10)") " PAO|", "MACHINE_LEARNING%DESCRIPTOR", id2str(pao%ml_descriptor)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,A11)") " PAO|", "MACHINE_LEARNING%DESCRIPTOR", id2str(pao%ml_descriptor)
|
||||
|
||||
CALL section_vals_val_get(pao_section, "MACHINE_LEARNING%TOLERANCE", r_val=pao%ml_tolerance)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "MACHINE_LEARNING%TOLERANCE", pao%ml_tolerance
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "MACHINE_LEARNING%TOLERANCE", pao%ml_tolerance
|
||||
|
||||
CALL section_vals_val_get(pao_section, "MACHINE_LEARNING%GP_NOISE_VAR", r_val=pao%gp_noise_var)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "MACHINE_LEARNING%GP_NOISE_VAR", pao%gp_noise_var
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "MACHINE_LEARNING%GP_NOISE_VAR", pao%gp_noise_var
|
||||
|
||||
CALL section_vals_val_get(pao_section, "MACHINE_LEARNING%GP_SCALE", r_val=pao%gp_scale)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,E10.1)") " PAO|", "MACHINE_LEARNING%GP_SCALE", pao%gp_scale
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,E11.1)") " PAO|", "MACHINE_LEARNING%GP_SCALE", pao%gp_scale
|
||||
|
||||
! parse MACHINE_LEARNING%TRAINING_SET section
|
||||
training_set_section => section_vals_get_subs_vals(pao_section, "MACHINE_LEARNING%TRAINING_SET")
|
||||
|
|
@ -198,7 +199,7 @@ CONTAINS
|
|||
DO i = 1, ntrainfiles
|
||||
CALL section_vals_val_get(training_set_section, "_DEFAULT_KEYWORD_", &
|
||||
i_rep_val=i, c_val=pao%ml_training_set(i)%fn)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T71,A)") " PAO|", "MACHINE_LEARNING%TRAINING_SET", &
|
||||
IF (pao%iw > 0) WRITE (pao%iw, "(A,T40,A,T70,A)") " PAO|", "MACHINE_LEARNING%TRAINING_SET", &
|
||||
TRIM(pao%ml_training_set(i)%fn)
|
||||
END DO
|
||||
|
||||
|
|
@ -213,7 +214,7 @@ CONTAINS
|
|||
! **************************************************************************************************
|
||||
FUNCTION id2str(id) RESULT(s)
|
||||
INTEGER :: id
|
||||
CHARACTER(LEN=10) :: s
|
||||
CHARACTER(LEN=11) :: s
|
||||
|
||||
SELECT CASE (id)
|
||||
CASE (pao_gth_param)
|
||||
|
|
@ -224,6 +225,8 @@ CONTAINS
|
|||
s = "FOCK"
|
||||
CASE (pao_exp_param)
|
||||
s = "EXP"
|
||||
CASE (pao_equi_param)
|
||||
s = "EQUIVARIANT"
|
||||
CASE (pao_opt_cg)
|
||||
s = "CG"
|
||||
CASE (pao_opt_bfgs)
|
||||
|
|
@ -296,12 +299,13 @@ CONTAINS
|
|||
! Parametrization **********************************************************
|
||||
CALL keyword_create(keyword, __LOCATION__, name="PARAMETERIZATION", &
|
||||
description="Parametrization of the mapping between the primary and the PAO basis.", &
|
||||
enum_c_vals=s2a("ROTINV", "FOCK", "GTH", "EXP"), &
|
||||
enum_i_vals=(/pao_rotinv_param, pao_fock_param, pao_gth_param, pao_exp_param/), &
|
||||
enum_c_vals=s2a("ROTINV", "FOCK", "GTH", "EXP", "EQUIVARIANT"), &
|
||||
enum_i_vals=(/pao_rotinv_param, pao_fock_param, pao_gth_param, pao_exp_param, pao_equi_param/), &
|
||||
enum_desc=s2a("Rotational invariant parametrization (machine learnable)", &
|
||||
"Fock matrix parametrization", &
|
||||
"Parametrization based on GTH pseudo potentials", &
|
||||
"Original matrix exponential parametrization"), &
|
||||
"Original matrix exponential parametrization", &
|
||||
"Equivariant parametrization"), &
|
||||
default_i_val=pao_rotinv_param)
|
||||
CALL section_add_keyword(section, keyword)
|
||||
CALL keyword_release(keyword)
|
||||
|
|
|
|||
|
|
@ -512,7 +512,7 @@ CONTAINS
|
|||
WRITE (unit_nr, "(A,5X,I0)") "Version", file_format_version
|
||||
WRITE (unit_nr, "(A,5X,F20.10)") "Energy", energy
|
||||
WRITE (unit_nr, "(A,5X,I0)") "Step", pao%istep
|
||||
WRITE (unit_nr, "(A,5X,A10)") "Parametrization", id2str(pao%parameterization)
|
||||
WRITE (unit_nr, "(A,5X,A)") "Parametrization", id2str(pao%parameterization)
|
||||
|
||||
! write kinds
|
||||
WRITE (unit_nr, "(A,5X,I0)") "Nkinds", SIZE(atomic_kind_set)
|
||||
|
|
|
|||
|
|
@ -39,18 +39,17 @@ MODULE pao_main
|
|||
USE pao_methods, ONLY: &
|
||||
pao_add_forces, pao_build_core_hamiltonian, pao_build_diag_distribution, &
|
||||
pao_build_matrix_X, pao_build_orthogonalizer, pao_build_selector, pao_calc_energy, &
|
||||
pao_calc_outer_grad_lnv, pao_check_grad, pao_check_trace_ps, pao_guess_initial_P, &
|
||||
pao_init_kinds, pao_print_atom_info, pao_store_P, pao_test_convergence
|
||||
pao_check_grad, pao_check_trace_ps, pao_guess_initial_P, pao_init_kinds, &
|
||||
pao_print_atom_info, pao_store_P, pao_test_convergence
|
||||
USE pao_ml, ONLY: pao_ml_init,&
|
||||
pao_ml_predict
|
||||
USE pao_optimizer, ONLY: pao_opt_finalize,&
|
||||
pao_opt_init,&
|
||||
pao_opt_new_dir
|
||||
USE pao_param, ONLY: pao_calc_U,&
|
||||
USE pao_param, ONLY: pao_calc_AB,&
|
||||
pao_param_finalize,&
|
||||
pao_param_init,&
|
||||
pao_param_initial_guess,&
|
||||
pao_update_AB
|
||||
pao_param_initial_guess
|
||||
USE pao_types, ONLY: pao_env_type
|
||||
USE qs_environment_types, ONLY: get_qs_env,&
|
||||
qs_environment_type
|
||||
|
|
@ -112,7 +111,6 @@ CONTAINS
|
|||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_optimization_start'
|
||||
|
||||
INTEGER :: handle
|
||||
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
|
||||
TYPE(ls_mstruct_type), POINTER :: ls_mstruct
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(section_vals_type), POINTER :: input, section
|
||||
|
|
@ -120,10 +118,7 @@ CONTAINS
|
|||
IF (.NOT. ls_scf_env%do_pao) RETURN
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
CALL get_qs_env(qs_env=qs_env, &
|
||||
matrix_s=matrix_s, &
|
||||
input=input)
|
||||
|
||||
CALL get_qs_env(qs_env, input=input)
|
||||
pao => ls_scf_env%pao_env
|
||||
ls_mstruct => ls_scf_env%ls_mstruct
|
||||
|
||||
|
|
@ -169,20 +164,13 @@ CONTAINS
|
|||
CALL dbcsr_copy(pao%matrix_G, pao%matrix_X)
|
||||
CALL dbcsr_set(pao%matrix_G, 0.0_dp)
|
||||
|
||||
CALL dbcsr_create(pao%matrix_U, &
|
||||
name="PAO matrix_U", &
|
||||
matrix_type="N", &
|
||||
dist=pao%diag_distribution, &
|
||||
template=matrix_s(1)%matrix)
|
||||
CALL dbcsr_reserve_diag_blocks(pao%matrix_U)
|
||||
|
||||
CALL dbcsr_create(ls_mstruct%matrix_A, template=pao%matrix_Y)
|
||||
CALL dbcsr_reserve_diag_blocks(ls_mstruct%matrix_A)
|
||||
CALL dbcsr_create(ls_mstruct%matrix_B, template=pao%matrix_Y)
|
||||
CALL dbcsr_reserve_diag_blocks(ls_mstruct%matrix_B)
|
||||
|
||||
! fill PAO transformation matrices
|
||||
CALL pao_update_AB(pao, qs_env, ls_mstruct)
|
||||
CALL pao_calc_AB(pao, qs_env, ls_scf_env, gradient=.FALSE.)
|
||||
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_optimization_start
|
||||
|
|
@ -203,7 +191,7 @@ CONTAINS
|
|||
INTEGER :: handle, icycle
|
||||
LOGICAL :: cycle_converged, do_mixing, should_stop
|
||||
REAL(KIND=dp) :: energy, penalty
|
||||
TYPE(dbcsr_type) :: matrix_M, matrix_X_mixing
|
||||
TYPE(dbcsr_type) :: matrix_X_mixing
|
||||
TYPE(ls_mstruct_type), POINTER :: ls_mstruct
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
|
||||
|
|
@ -258,9 +246,7 @@ CONTAINS
|
|||
IF (pao%linesearch%starts) THEN
|
||||
icycle = icycle + 1
|
||||
! calc new gradient including penalty terms
|
||||
CALL pao_calc_outer_grad_lnv(qs_env, ls_scf_env, matrix_M)
|
||||
CALL pao_calc_U(pao, qs_env, matrix_M=matrix_M, matrix_G=pao%matrix_G, penalty=penalty)
|
||||
CALL dbcsr_release(matrix_M)
|
||||
CALL pao_calc_AB(pao, qs_env, ls_scf_env, gradient=.TRUE., penalty=penalty)
|
||||
CALL pao_check_grad(pao, qs_env, ls_scf_env)
|
||||
|
||||
! calculate new direction for line-search
|
||||
|
|
@ -376,7 +362,6 @@ CONTAINS
|
|||
! We keep pao%matrix_X for next scf-run, e.g. during MD or GEO-OPT
|
||||
CALL dbcsr_release(pao%matrix_X_orig)
|
||||
CALL dbcsr_release(pao%matrix_G)
|
||||
CALL dbcsr_release(pao%matrix_U)
|
||||
CALL dbcsr_release(ls_mstruct%matrix_A)
|
||||
CALL dbcsr_release(ls_mstruct%matrix_B)
|
||||
|
||||
|
|
|
|||
|
|
@ -27,13 +27,12 @@ MODULE pao_methods
|
|||
dbcsr_copy_into_existing, dbcsr_create, dbcsr_desymmetrize, dbcsr_distribution_get, &
|
||||
dbcsr_distribution_new, dbcsr_distribution_type, dbcsr_dot, dbcsr_filter, &
|
||||
dbcsr_get_block_p, dbcsr_get_info, dbcsr_iterator_blocks_left, dbcsr_iterator_next_block, &
|
||||
dbcsr_iterator_start, dbcsr_iterator_stop, dbcsr_iterator_type, dbcsr_multiply, &
|
||||
dbcsr_p_type, dbcsr_release, dbcsr_reserve_diag_blocks, dbcsr_scale, dbcsr_set, dbcsr_type
|
||||
dbcsr_iterator_start, dbcsr_iterator_stop, dbcsr_iterator_type, dbcsr_p_type, &
|
||||
dbcsr_release, dbcsr_reserve_diag_blocks, dbcsr_scale, dbcsr_set, dbcsr_type
|
||||
USE dm_ls_scf_methods, ONLY: density_matrix_trs4,&
|
||||
ls_scf_init_matrix_S
|
||||
USE dm_ls_scf_qs, ONLY: ls_scf_dm_to_ks,&
|
||||
ls_scf_qs_atomic_guess,&
|
||||
matrix_decluster,&
|
||||
matrix_ls_to_qs,&
|
||||
matrix_qs_to_ls
|
||||
USE dm_ls_scf_types, ONLY: ls_mstruct_type,&
|
||||
|
|
@ -47,9 +46,8 @@ MODULE pao_methods
|
|||
USE message_passing, ONLY: mp_max,&
|
||||
mp_sum
|
||||
USE pao_ml, ONLY: pao_ml_forces
|
||||
USE pao_param, ONLY: pao_calc_U,&
|
||||
pao_param_count,&
|
||||
pao_update_AB
|
||||
USE pao_param, ONLY: pao_calc_AB,&
|
||||
pao_param_count
|
||||
USE pao_types, ONLY: pao_env_type
|
||||
USE particle_types, ONLY: particle_type
|
||||
USE qs_energy_types, ONLY: qs_energy_type
|
||||
|
|
@ -83,7 +81,6 @@ MODULE pao_methods
|
|||
PUBLIC :: pao_test_convergence
|
||||
PUBLIC :: pao_calc_energy, pao_check_trace_ps
|
||||
PUBLIC :: pao_store_P, pao_add_forces, pao_guess_initial_P
|
||||
PUBLIC :: pao_calc_outer_grad_lnv
|
||||
PUBLIC :: pao_check_grad
|
||||
|
||||
CONTAINS
|
||||
|
|
@ -228,13 +225,19 @@ CONTAINS
|
|||
CALL dbcsr_iterator_stop(iter)
|
||||
!$OMP END PARALLEL
|
||||
|
||||
! store a copy that is distributed according to pao%diag_distribution
|
||||
! store a copies of N and N_inv that are distributed according to pao%diag_distribution
|
||||
CALL dbcsr_create(pao%matrix_N_diag, &
|
||||
name="PAO matrix_N_diag", &
|
||||
dist=pao%diag_distribution, &
|
||||
template=matrix_s(1)%matrix)
|
||||
CALL dbcsr_reserve_diag_blocks(pao%matrix_N_diag)
|
||||
CALL dbcsr_complete_redistribute(pao%matrix_N, pao%matrix_N_diag)
|
||||
CALL dbcsr_create(pao%matrix_N_inv_diag, &
|
||||
name="PAO matrix_N_inv_diag", &
|
||||
dist=pao%diag_distribution, &
|
||||
template=matrix_s(1)%matrix)
|
||||
CALL dbcsr_reserve_diag_blocks(pao%matrix_N_inv_diag)
|
||||
CALL dbcsr_complete_redistribute(pao%matrix_N_inv, pao%matrix_N_inv_diag)
|
||||
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_build_orthogonalizer
|
||||
|
|
@ -491,7 +494,7 @@ CONTAINS
|
|||
CALL timeset(routineN, handle)
|
||||
|
||||
! calculate matrix U, which determines the pao basis
|
||||
CALL pao_update_AB(pao, qs_env, ls_scf_env%ls_mstruct, penalty=penalty)
|
||||
CALL pao_calc_AB(pao, qs_env, ls_scf_env, gradient=.FALSE., penalty=penalty)
|
||||
|
||||
! calculat S, S_inv, S_sqrt, and S_sqrt_inv in the new pao basis
|
||||
CALL pao_rebuild_S(qs_env, ls_scf_env)
|
||||
|
|
@ -671,234 +674,6 @@ CONTAINS
|
|||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_dm_trs4
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Helper routine, calculates partial derivative dE/dU
|
||||
!> \param qs_env ...
|
||||
!> \param ls_scf_env ...
|
||||
!> \param matrix_M_diag the derivate, matrix uses pao%diag_distribution
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_outer_grad_lnv(qs_env, ls_scf_env, matrix_M_diag)
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
TYPE(ls_scf_env_type), TARGET :: ls_scf_env
|
||||
TYPE(dbcsr_type) :: matrix_M_diag
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_outer_grad_lnv'
|
||||
|
||||
INTEGER :: handle, nspin
|
||||
INTEGER, DIMENSION(:), POINTER :: pao_blk_sizes
|
||||
REAL(KIND=dp) :: filter_eps
|
||||
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s, rho_ao
|
||||
TYPE(dbcsr_type) :: matrix_HB, matrix_HPS, matrix_M, matrix_M1, matrix_M1_dc, matrix_M2, &
|
||||
matrix_M2_dc, matrix_M3, matrix_M3_dc, matrix_NHB, matrix_NHBM2, matrix_NPA, &
|
||||
matrix_NPAM1, matrix_NSB, matrix_NSBM3, matrix_PA, matrix_PH, matrix_PHP, matrix_PSP, &
|
||||
matrix_SB, matrix_SP
|
||||
TYPE(dft_control_type), POINTER :: dft_control
|
||||
TYPE(ls_mstruct_type), POINTER :: ls_mstruct
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_rho_type), POINTER :: rho
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
|
||||
ls_mstruct => ls_scf_env%ls_mstruct
|
||||
pao => ls_scf_env%pao_env
|
||||
|
||||
CALL get_qs_env(qs_env, &
|
||||
rho=rho, &
|
||||
matrix_ks=matrix_ks, &
|
||||
matrix_s=matrix_s, &
|
||||
dft_control=dft_control)
|
||||
CALL qs_rho_get(rho, rho_ao=rho_ao)
|
||||
nspin = dft_control%nspins
|
||||
filter_eps = ls_scf_env%eps_filter
|
||||
|
||||
CALL dbcsr_get_info(ls_mstruct%matrix_A, col_blk_size=pao_blk_sizes)
|
||||
|
||||
IF (nspin /= 1) CPABORT("open shell not yet implemented")
|
||||
!TODO: handle openshell case properly
|
||||
|
||||
! notation according to pao_math_lnv.pdf
|
||||
|
||||
! calculation uses distribution of matrix_s, after we redistribute using pao%diag_distribution
|
||||
CALL dbcsr_create(matrix_M, template=matrix_s(1)%matrix, matrix_type="N")
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_M)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! calculate need products in pao basis
|
||||
CALL dbcsr_create(matrix_PH, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, ls_scf_env%matrix_p(1), ls_scf_env%matrix_ks(1), &
|
||||
0.0_dp, matrix_PH, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_PHP, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_PH, ls_scf_env%matrix_p(1), &
|
||||
0.0_dp, matrix_PHP, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_SP, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, ls_scf_env%matrix_s, ls_scf_env%matrix_p(1), &
|
||||
0.0_dp, matrix_SP, filter_eps=filter_eps)
|
||||
|
||||
IF (SIZE(ls_scf_env%matrix_p) == 1) CALL dbcsr_scale(matrix_SP, 0.5_dp)
|
||||
|
||||
CALL dbcsr_create(matrix_HPS, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "T", 1.0_dp, ls_scf_env%matrix_ks(1), matrix_SP, &
|
||||
0.0_dp, matrix_HPS, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_PSP, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, ls_scf_env%matrix_p(1), matrix_SP, &
|
||||
0.0_dp, matrix_PSP, filter_eps=filter_eps)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! M1 = dE_lnv / dP_pao
|
||||
CALL dbcsr_create(matrix_M1, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
|
||||
CALL dbcsr_multiply("N", "T", 3.0_dp, ls_scf_env%matrix_ks(1), matrix_SP, &
|
||||
1.0_dp, matrix_M1, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "N", 3.0_dp, matrix_SP, ls_scf_env%matrix_ks(1), &
|
||||
1.0_dp, matrix_M1, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "T", -2.0_dp, matrix_HPS, matrix_SP, &
|
||||
1.0_dp, matrix_M1, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "N", -2.0_dp, matrix_SP, matrix_HPS, &
|
||||
1.0_dp, matrix_M1, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "T", -2.0_dp, matrix_SP, matrix_HPS, &
|
||||
1.0_dp, matrix_M1, filter_eps=filter_eps)
|
||||
|
||||
! reverse possible molecular clustering
|
||||
CALL dbcsr_create(matrix_M1_dc, &
|
||||
template=matrix_s(1)%matrix, &
|
||||
row_blk_size=pao_blk_sizes, &
|
||||
col_blk_size=pao_blk_sizes)
|
||||
CALL matrix_decluster(matrix_M1_dc, matrix_M1, ls_mstruct)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! M2 = dE_lnv / dH
|
||||
CALL dbcsr_create(matrix_M2, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
|
||||
CALL dbcsr_add(matrix_M2, matrix_PSP, 1.0_dp, 3.0_dp)
|
||||
|
||||
CALL dbcsr_multiply("N", "N", -2.0_dp, matrix_PSP, matrix_SP, &
|
||||
1.0_dp, matrix_M2, filter_eps=filter_eps)
|
||||
|
||||
! reverse possible molecular clustering
|
||||
CALL dbcsr_create(matrix_M2_dc, &
|
||||
template=matrix_s(1)%matrix, &
|
||||
row_blk_size=pao_blk_sizes, &
|
||||
col_blk_size=pao_blk_sizes)
|
||||
CALL matrix_decluster(matrix_M2_dc, matrix_M2, ls_mstruct)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! M3 = dE_lnv / dS
|
||||
CALL dbcsr_create(matrix_M3, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
|
||||
CALL dbcsr_add(matrix_M3, matrix_PHP, 1.0_dp, 3.0_dp)
|
||||
|
||||
CALL dbcsr_multiply("N", "N", -2.0_dp, matrix_PHP, matrix_SP, &
|
||||
1.0_dp, matrix_M3, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "T", -2.0_dp, matrix_PSP, matrix_PH, &
|
||||
1.0_dp, matrix_M3, filter_eps=filter_eps)
|
||||
|
||||
! reverse possible molecular clustering
|
||||
CALL dbcsr_create(matrix_M3_dc, &
|
||||
template=matrix_s(1)%matrix, &
|
||||
row_blk_size=pao_blk_sizes, &
|
||||
col_blk_size=pao_blk_sizes)
|
||||
CALL matrix_decluster(matrix_M3_dc, matrix_M3, ls_mstruct)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! combine M1 with matrices from primary basis
|
||||
CALL dbcsr_create(matrix_PA, template=ls_mstruct%matrix_A, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, rho_ao(1)%matrix, ls_mstruct%matrix_A, &
|
||||
0.0_dp, matrix_PA, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_NPA, template=ls_mstruct%matrix_A, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, pao%matrix_N_inv, matrix_PA, &
|
||||
0.0_dp, matrix_NPA, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_NPAM1, template=ls_mstruct%matrix_A, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_NPA, matrix_M1_dc, &
|
||||
0.0_dp, matrix_NPAM1, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "T", 1.0_dp, matrix_NPAM1, pao%matrix_Y, &
|
||||
1.0_dp, matrix_M, filter_eps=filter_eps)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! combine M2 with matrices from primary basis
|
||||
CALL dbcsr_create(matrix_HB, template=ls_mstruct%matrix_B, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_ks(1)%matrix, ls_mstruct%matrix_B, &
|
||||
0.0_dp, matrix_HB, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_NHB, template=ls_mstruct%matrix_B, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, pao%matrix_N, matrix_HB, &
|
||||
0.0_dp, matrix_NHB, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_NHBM2, template=ls_mstruct%matrix_B, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_NHB, matrix_M2_dc, &
|
||||
0.0_dp, matrix_NHBM2, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "T", 1.0_dp, matrix_NHBM2, pao%matrix_Y, &
|
||||
1.0_dp, matrix_M, filter_eps=filter_eps)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! combine M3 with matrices from primary basis
|
||||
CALL dbcsr_create(matrix_SB, template=ls_mstruct%matrix_B, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_s(1)%matrix, ls_mstruct%matrix_B, &
|
||||
0.0_dp, matrix_SB, filter_eps=filter_eps)
|
||||
|
||||
IF (SIZE(ls_scf_env%matrix_p) == 1) CALL dbcsr_scale(matrix_SB, 0.5_dp)
|
||||
|
||||
CALL dbcsr_create(matrix_NSB, template=ls_mstruct%matrix_B, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, pao%matrix_N, matrix_SB, &
|
||||
0.0_dp, matrix_NSB, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_NSBM3, template=ls_mstruct%matrix_B, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_NSB, matrix_M3_dc, &
|
||||
0.0_dp, matrix_NSBM3, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "T", 1.0_dp, matrix_NSBM3, pao%matrix_Y, &
|
||||
1.0_dp, matrix_M, filter_eps=filter_eps)
|
||||
|
||||
IF (SIZE(ls_scf_env%matrix_p) == 1) CALL dbcsr_scale(matrix_M, 2.0_dp)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! redistribute using pao%diag_distribution
|
||||
CALL dbcsr_create(matrix_M_diag, &
|
||||
name="PAO matrix_M", &
|
||||
matrix_type="N", &
|
||||
dist=pao%diag_distribution, &
|
||||
template=matrix_s(1)%matrix)
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_M_diag)
|
||||
CALL dbcsr_complete_redistribute(matrix_M, matrix_M_diag)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! cleanup: TODO release matrices as early as possible
|
||||
CALL dbcsr_release(matrix_PH)
|
||||
CALL dbcsr_release(matrix_PHP)
|
||||
CALL dbcsr_release(matrix_SP)
|
||||
CALL dbcsr_release(matrix_HPS)
|
||||
CALL dbcsr_release(matrix_PSP)
|
||||
CALL dbcsr_release(matrix_M)
|
||||
CALL dbcsr_release(matrix_M1)
|
||||
CALL dbcsr_release(matrix_M2)
|
||||
CALL dbcsr_release(matrix_M3)
|
||||
CALL dbcsr_release(matrix_M1_dc)
|
||||
CALL dbcsr_release(matrix_M2_dc)
|
||||
CALL dbcsr_release(matrix_M3_dc)
|
||||
CALL dbcsr_release(matrix_PA)
|
||||
CALL dbcsr_release(matrix_NPA)
|
||||
CALL dbcsr_release(matrix_NPAM1)
|
||||
CALL dbcsr_release(matrix_HB)
|
||||
CALL dbcsr_release(matrix_NHB)
|
||||
CALL dbcsr_release(matrix_NHBM2)
|
||||
CALL dbcsr_release(matrix_SB)
|
||||
CALL dbcsr_release(matrix_NSB)
|
||||
CALL dbcsr_release(matrix_NSBM3)
|
||||
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_calc_outer_grad_lnv
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Debugging routine for checking the analytic gradient.
|
||||
!> \param pao ...
|
||||
|
|
@ -1189,7 +964,6 @@ CONTAINS
|
|||
INTEGER :: handle, iatom, natoms
|
||||
REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: forces
|
||||
TYPE(cp_para_env_type), POINTER :: para_env
|
||||
TYPE(dbcsr_type) :: matrix_M
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
|
||||
|
||||
|
|
@ -1213,10 +987,7 @@ CONTAINS
|
|||
natom=natoms)
|
||||
|
||||
ALLOCATE (forces(natoms, 3))
|
||||
forces(:, :) = 0.0_dp
|
||||
CALL pao_calc_outer_grad_lnv(qs_env, ls_scf_env, matrix_M)
|
||||
CALL pao_calc_U(pao, qs_env, matrix_M, pao%matrix_G, forces=forces) ! without penalty terms
|
||||
CALL dbcsr_release(matrix_M)
|
||||
CALL pao_calc_AB(pao, qs_env, ls_scf_env, gradient=.TRUE., forces=forces) ! without penalty terms
|
||||
|
||||
IF (SIZE(pao%ml_training_set) > 0) &
|
||||
CALL pao_ml_forces(pao, qs_env, pao%matrix_G, forces)
|
||||
|
|
|
|||
215
src/pao_param.F
215
src/pao_param.F
|
|
@ -10,38 +10,40 @@
|
|||
!> \author Ole Schuett
|
||||
! **************************************************************************************************
|
||||
MODULE pao_param
|
||||
USE cp_log_handling, ONLY: cp_to_string
|
||||
USE dbcsr_api, ONLY: &
|
||||
dbcsr_complete_redistribute, dbcsr_copy, dbcsr_create, dbcsr_frobenius_norm, &
|
||||
dbcsr_get_block_p, dbcsr_get_info, dbcsr_iterator_blocks_left, dbcsr_iterator_next_block, &
|
||||
dbcsr_iterator_start, dbcsr_iterator_stop, dbcsr_iterator_type, dbcsr_multiply, &
|
||||
dbcsr_p_type, dbcsr_release, dbcsr_reserve_diag_blocks, dbcsr_type
|
||||
USE dm_ls_scf_types, ONLY: ls_mstruct_type
|
||||
USE dbcsr_api, ONLY: dbcsr_copy,&
|
||||
dbcsr_frobenius_norm,&
|
||||
dbcsr_multiply,&
|
||||
dbcsr_release,&
|
||||
dbcsr_type
|
||||
USE dm_ls_scf_types, ONLY: ls_scf_env_type
|
||||
USE kinds, ONLY: dp
|
||||
USE message_passing, ONLY: mp_comm_type,&
|
||||
mp_max
|
||||
USE pao_input, ONLY: pao_exp_param,&
|
||||
USE pao_input, ONLY: pao_equi_param,&
|
||||
pao_exp_param,&
|
||||
pao_fock_param,&
|
||||
pao_gth_param,&
|
||||
pao_rotinv_param
|
||||
USE pao_param_exp, ONLY: pao_calc_U_exp,&
|
||||
USE pao_param_equi, ONLY: pao_calc_AB_equi,&
|
||||
pao_param_count_equi,&
|
||||
pao_param_finalize_equi,&
|
||||
pao_param_init_equi,&
|
||||
pao_param_initguess_equi
|
||||
USE pao_param_exp, ONLY: pao_calc_AB_exp,&
|
||||
pao_param_count_exp,&
|
||||
pao_param_finalize_exp,&
|
||||
pao_param_init_exp,&
|
||||
pao_param_initguess_exp
|
||||
USE pao_param_gth, ONLY: pao_calc_U_gth,&
|
||||
USE pao_param_gth, ONLY: pao_calc_AB_gth,&
|
||||
pao_param_count_gth,&
|
||||
pao_param_finalize_gth,&
|
||||
pao_param_init_gth,&
|
||||
pao_param_initguess_gth
|
||||
USE pao_param_linpot, ONLY: pao_calc_U_linpot,&
|
||||
USE pao_param_linpot, ONLY: pao_calc_AB_linpot,&
|
||||
pao_param_count_linpot,&
|
||||
pao_param_finalize_linpot,&
|
||||
pao_param_init_linpot,&
|
||||
pao_param_initguess_linpot
|
||||
USE pao_types, ONLY: pao_env_type
|
||||
USE qs_environment_types, ONLY: get_qs_env,&
|
||||
qs_environment_type
|
||||
USE qs_environment_types, ONLY: qs_environment_type
|
||||
#include "./base/base_uses.f90"
|
||||
|
||||
IMPLICIT NONE
|
||||
|
|
@ -50,74 +52,53 @@ MODULE pao_param
|
|||
|
||||
CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'pao_param'
|
||||
|
||||
PUBLIC :: pao_update_AB, pao_param_count, pao_param_initial_guess
|
||||
PUBLIC :: pao_param_init, pao_param_finalize, pao_calc_U
|
||||
PUBLIC :: pao_calc_AB, pao_param_count, pao_param_initial_guess
|
||||
PUBLIC :: pao_param_init, pao_param_finalize
|
||||
|
||||
CONTAINS
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Takes current matrix_X and recalculates derived matrices U, A, and B.
|
||||
!> \brief Takes current matrix_X and calculates the matrices A and B.
|
||||
!> \param pao ...
|
||||
!> \param qs_env ...
|
||||
!> \param ls_mstruct ...
|
||||
!> \param ls_scf_env ...
|
||||
!> \param gradient ...
|
||||
!> \param penalty ...
|
||||
!> \param forces ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_update_AB(pao, qs_env, ls_mstruct, penalty)
|
||||
SUBROUTINE pao_calc_AB(pao, qs_env, ls_scf_env, gradient, penalty, forces)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
TYPE(ls_mstruct_type) :: ls_mstruct
|
||||
TYPE(ls_scf_env_type), TARGET :: ls_scf_env
|
||||
LOGICAL, INTENT(IN) :: gradient
|
||||
REAL(dp), INTENT(OUT), OPTIONAL :: penalty
|
||||
REAL(dp), DIMENSION(:, :), INTENT(OUT), OPTIONAL :: forces
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_update_AB'
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_AB'
|
||||
|
||||
INTEGER :: acol, arow, handle, iatom
|
||||
LOGICAL :: found
|
||||
REAL(dp), DIMENSION(:, :), POINTER :: block_A, block_B, block_N, block_N_inv, &
|
||||
block_U, block_Y
|
||||
TYPE(dbcsr_iterator_type) :: iter
|
||||
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
|
||||
TYPE(dbcsr_type) :: matrix_U
|
||||
INTEGER :: handle
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
|
||||
CALL pao_calc_U(pao, qs_env, penalty=penalty) !update matrix_U = Function of matrix_X
|
||||
IF (PRESENT(penalty)) penalty = 0.0_dp
|
||||
IF (PRESENT(forces)) forces(:, :) = 0.0_dp
|
||||
|
||||
! pao%matrix_U uses pao%diag_distribution, need to redistribute using distribution of matrix_s
|
||||
CALL get_qs_env(qs_env, matrix_s=matrix_s)
|
||||
CALL dbcsr_create(matrix_U, matrix_type="N", template=matrix_s(1)%matrix)
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_U)
|
||||
CALL dbcsr_complete_redistribute(pao%matrix_U, matrix_U)
|
||||
|
||||
! Multiplying diagonal matrices is a local operation.
|
||||
! To take advantage of this we're using an iterator instead of calling dbcsr_multiply().
|
||||
!$OMP PARALLEL DEFAULT(NONE) SHARED(pao,ls_mstruct,matrix_U) &
|
||||
!$OMP PRIVATE(iter,arow,acol,iatom,block_U,block_Y,block_A,block_B,block_N,block_N_inv,found)
|
||||
CALL dbcsr_iterator_start(iter, matrix_U)
|
||||
DO WHILE (dbcsr_iterator_blocks_left(iter))
|
||||
CALL dbcsr_iterator_next_block(iter, arow, acol, block_U)
|
||||
iatom = arow; CPASSERT(arow == acol)
|
||||
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_Y, row=iatom, col=iatom, block=block_Y, found=found)
|
||||
CPASSERT(ASSOCIATED(block_Y))
|
||||
|
||||
CALL dbcsr_get_block_p(matrix=ls_mstruct%matrix_A, row=iatom, col=iatom, block=block_A, found=found)
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_N_inv, row=iatom, col=iatom, block=block_N_inv, found=found)
|
||||
CPASSERT(ASSOCIATED(block_A) .AND. ASSOCIATED(block_N_inv))
|
||||
|
||||
CALL dbcsr_get_block_p(matrix=ls_mstruct%matrix_B, row=iatom, col=iatom, block=block_B, found=found)
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_N, row=iatom, col=iatom, block=block_N, found=found)
|
||||
CPASSERT(ASSOCIATED(block_B) .AND. ASSOCIATED(block_N))
|
||||
|
||||
block_A = MATMUL(MATMUL(block_N_inv, block_U), block_Y)
|
||||
block_B = MATMUL(MATMUL(block_N, block_U), block_Y)
|
||||
END DO
|
||||
CALL dbcsr_iterator_stop(iter)
|
||||
!$OMP END PARALLEL
|
||||
|
||||
CALL dbcsr_release(matrix_U)
|
||||
!calculate matrix_A/B = Function of matrix_X
|
||||
SELECT CASE (pao%parameterization)
|
||||
CASE (pao_exp_param)
|
||||
CALL pao_calc_AB_exp(pao, qs_env, ls_scf_env, gradient)
|
||||
CASE (pao_fock_param, pao_rotinv_param)
|
||||
CALL pao_calc_AB_linpot(pao, qs_env, ls_scf_env, gradient, penalty, forces)
|
||||
CASE (pao_gth_param)
|
||||
CALL pao_calc_AB_gth(pao, qs_env, ls_scf_env, gradient, penalty)
|
||||
CASE (pao_equi_param)
|
||||
CALL pao_calc_AB_equi(pao, qs_env, ls_scf_env, gradient, penalty)
|
||||
CASE DEFAULT
|
||||
CPABORT("PAO: unkown parametrization")
|
||||
END SELECT
|
||||
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_update_AB
|
||||
END SUBROUTINE pao_calc_AB
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Initialize PAO parametrization
|
||||
|
|
@ -141,6 +122,8 @@ CONTAINS
|
|||
CALL pao_param_init_linpot(pao, qs_env)
|
||||
CASE (pao_gth_param)
|
||||
CALL pao_param_init_gth(pao, qs_env)
|
||||
CASE (pao_equi_param)
|
||||
CALL pao_param_init_equi(pao)
|
||||
CASE DEFAULT
|
||||
CPABORT("PAO: unknown parametrization")
|
||||
END SELECT
|
||||
|
|
@ -169,6 +152,8 @@ CONTAINS
|
|||
CALL pao_param_finalize_linpot(pao)
|
||||
CASE (pao_gth_param)
|
||||
CALL pao_param_finalize_gth(pao)
|
||||
CASE (pao_equi_param)
|
||||
CALL pao_param_finalize_equi()
|
||||
CASE DEFAULT
|
||||
CPABORT("PAO: unknown parametrization")
|
||||
END SELECT
|
||||
|
|
@ -203,6 +188,8 @@ CONTAINS
|
|||
CALL pao_param_count_linpot(pao, qs_env, ikind=ikind, nparams=nparams)
|
||||
CASE (pao_gth_param)
|
||||
CALL pao_param_count_gth(qs_env, ikind=ikind, nparams=nparams)
|
||||
CASE (pao_equi_param)
|
||||
CALL pao_param_count_equi(qs_env, ikind=ikind, nparams=nparams)
|
||||
CASE DEFAULT
|
||||
CPABORT("PAO: unknown parametrization")
|
||||
END SELECT
|
||||
|
|
@ -235,6 +222,8 @@ CONTAINS
|
|||
CALL pao_param_initguess_linpot(pao, qs_env)
|
||||
CASE (pao_gth_param)
|
||||
CALL pao_param_initguess_gth(pao)
|
||||
CASE (pao_equi_param)
|
||||
CALL pao_param_initguess_equi(pao, qs_env)
|
||||
CASE DEFAULT
|
||||
CPABORT("PAO: unknown parametrization")
|
||||
END SELECT
|
||||
|
|
@ -254,102 +243,4 @@ CONTAINS
|
|||
|
||||
END SUBROUTINE pao_param_initial_guess
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Calculate new matrix U
|
||||
!> \param pao ...
|
||||
!> \param qs_env ...
|
||||
!> \param matrix_M ...
|
||||
!> \param matrix_G ...
|
||||
!> \param penalty ...
|
||||
!> \param forces ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_U(pao, qs_env, matrix_M, matrix_G, penalty, forces)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
TYPE(dbcsr_type), OPTIONAL :: matrix_M, matrix_G
|
||||
REAL(dp), INTENT(OUT), OPTIONAL :: penalty
|
||||
REAL(dp), DIMENSION(:, :), INTENT(INOUT), OPTIONAL :: forces
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_U'
|
||||
|
||||
INTEGER :: handle
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
|
||||
IF (PRESENT(penalty)) penalty = 0.0_dp
|
||||
|
||||
SELECT CASE (pao%parameterization)
|
||||
CASE (pao_exp_param)
|
||||
CALL pao_calc_U_exp(pao, matrix_M, matrix_G)
|
||||
CASE (pao_fock_param, pao_rotinv_param)
|
||||
CALL pao_calc_U_linpot(pao, qs_env, penalty, matrix_M, matrix_G, forces)
|
||||
CASE (pao_gth_param)
|
||||
CALL pao_calc_U_gth(pao, penalty, matrix_M, matrix_G)
|
||||
CASE DEFAULT
|
||||
CPABORT("PAO: unknown parametrization")
|
||||
END SELECT
|
||||
|
||||
CALL pao_assert_unitary(pao)
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_calc_U
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Debugging routine, check unitaryness of U
|
||||
!> \param pao ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_assert_unitary(pao)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_assert_unitary'
|
||||
|
||||
INTEGER :: acol, arow, group_handle, handle, i, &
|
||||
iatom, M, N
|
||||
INTEGER, DIMENSION(:), POINTER :: blk_sizes_pao, blk_sizes_pri
|
||||
REAL(dp) :: delta_max
|
||||
REAL(dp), DIMENSION(:, :), POINTER :: block_test, tmp1, tmp2
|
||||
TYPE(dbcsr_iterator_type) :: iter
|
||||
TYPE(mp_comm_type) :: group
|
||||
|
||||
IF (pao%check_unitary_tol < 0.0_dp) RETURN ! no checking
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
delta_max = 0.0_dp
|
||||
|
||||
CALL dbcsr_get_info(pao%matrix_Y, row_blk_size=blk_sizes_pri, col_blk_size=blk_sizes_pao)
|
||||
|
||||
!$OMP PARALLEL DEFAULT(NONE) SHARED(pao,blk_sizes_pri,blk_sizes_pao,delta_max) &
|
||||
!$OMP PRIVATE(iter,arow,acol,iatom,N,M,block_test,tmp1,tmp2)
|
||||
CALL dbcsr_iterator_start(iter, pao%matrix_U)
|
||||
DO WHILE (dbcsr_iterator_blocks_left(iter))
|
||||
CALL dbcsr_iterator_next_block(iter, arow, acol, block_test)
|
||||
iatom = arow; CPASSERT(arow == acol)
|
||||
N = blk_sizes_pri(iatom) ! size of primary basis
|
||||
M = blk_sizes_pao(iatom) ! size of pao basis
|
||||
|
||||
! we only need the upper left "PAO-corner" to be unitary
|
||||
ALLOCATE (tmp1(N, M), tmp2(M, M))
|
||||
tmp1 = block_test(:, 1:M)
|
||||
tmp2 = MATMUL(TRANSPOSE(tmp1), tmp1)
|
||||
DO i = 1, M
|
||||
tmp2(i, i) = tmp2(i, i) - 1.0_dp
|
||||
END DO
|
||||
|
||||
!$OMP ATOMIC
|
||||
delta_max = MAX(delta_max, MAXVAL(ABS(tmp2)))
|
||||
|
||||
DEALLOCATE (tmp1, tmp2)
|
||||
END DO
|
||||
CALL dbcsr_iterator_stop(iter)
|
||||
!$OMP END PARALLEL
|
||||
|
||||
CALL dbcsr_get_info(pao%matrix_U, group=group_handle)
|
||||
CALL group%set_handle(group_handle)
|
||||
CALL mp_max(delta_max, group)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, *) 'PAO| checked unitaryness, max delta:', delta_max
|
||||
IF (delta_max > pao%check_unitary_tol) &
|
||||
CPABORT("Found bad unitaryness:"//cp_to_string(delta_max))
|
||||
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_assert_unitary
|
||||
|
||||
END MODULE pao_param
|
||||
|
|
|
|||
317
src/pao_param_equi.F
Normal file
317
src/pao_param_equi.F
Normal file
|
|
@ -0,0 +1,317 @@
|
|||
!--------------------------------------------------------------------------------------------------!
|
||||
! CP2K: A general program to perform molecular dynamics simulations !
|
||||
! Copyright 2000-2023 CP2K developers group <https://cp2k.org> !
|
||||
! !
|
||||
! SPDX-License-Identifier: GPL-2.0-or-later !
|
||||
!--------------------------------------------------------------------------------------------------!
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Equivariant parametrization
|
||||
!> \author Ole Schuett
|
||||
! **************************************************************************************************
|
||||
MODULE pao_param_equi
|
||||
USE basis_set_types, ONLY: gto_basis_set_type
|
||||
USE dbcsr_api, ONLY: &
|
||||
dbcsr_get_block_p, dbcsr_get_info, dbcsr_iterator_blocks_left, dbcsr_iterator_next_block, &
|
||||
dbcsr_iterator_start, dbcsr_iterator_stop, dbcsr_iterator_type, dbcsr_release, dbcsr_type
|
||||
USE dm_ls_scf_types, ONLY: ls_mstruct_type,&
|
||||
ls_scf_env_type
|
||||
USE kinds, ONLY: dp
|
||||
USE mathlib, ONLY: diamat_all
|
||||
USE pao_param_methods, ONLY: pao_calc_grad_lnv_wrt_AB
|
||||
USE pao_potentials, ONLY: pao_guess_initial_potential
|
||||
USE pao_types, ONLY: pao_env_type
|
||||
USE qs_environment_types, ONLY: get_qs_env,&
|
||||
qs_environment_type
|
||||
USE qs_kind_types, ONLY: get_qs_kind,&
|
||||
qs_kind_type
|
||||
#include "./base/base_uses.f90"
|
||||
|
||||
IMPLICIT NONE
|
||||
|
||||
PRIVATE
|
||||
|
||||
CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'pao_param_equi'
|
||||
|
||||
PUBLIC :: pao_param_init_equi, pao_param_finalize_equi, pao_calc_AB_equi
|
||||
PUBLIC :: pao_param_count_equi, pao_param_initguess_equi
|
||||
|
||||
CONTAINS
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Initialize equivariant parametrization
|
||||
!> \param pao ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_param_init_equi(pao)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
|
||||
IF (pao%precondition) &
|
||||
CPABORT("PAO preconditioning not supported for selected parametrization.")
|
||||
|
||||
END SUBROUTINE pao_param_init_equi
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Finalize equivariant parametrization
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_param_finalize_equi()
|
||||
|
||||
! Nothing to do.
|
||||
|
||||
END SUBROUTINE pao_param_finalize_equi
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Returns the number of parameters for given atomic kind
|
||||
!> \param qs_env ...
|
||||
!> \param ikind ...
|
||||
!> \param nparams ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_param_count_equi(qs_env, ikind, nparams)
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
INTEGER, INTENT(IN) :: ikind
|
||||
INTEGER, INTENT(OUT) :: nparams
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_param_count_equi'
|
||||
|
||||
INTEGER :: pao_basis_size, pri_basis_size
|
||||
TYPE(gto_basis_set_type), POINTER :: basis_set
|
||||
TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
|
||||
|
||||
CALL get_qs_env(qs_env, qs_kind_set=qs_kind_set)
|
||||
CALL get_qs_kind(qs_kind_set(ikind), &
|
||||
basis_set=basis_set, &
|
||||
pao_basis_size=pao_basis_size)
|
||||
pri_basis_size = basis_set%nsgf
|
||||
|
||||
nparams = pao_basis_size*pri_basis_size
|
||||
|
||||
END SUBROUTINE pao_param_count_equi
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Fills matrix_X with an initial guess
|
||||
!> \param pao ...
|
||||
!> \param qs_env ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_param_initguess_equi(pao, qs_env)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_param_initguess_equi'
|
||||
|
||||
INTEGER :: acol, arow, handle, i, iatom, m, n
|
||||
INTEGER, DIMENSION(:), POINTER :: blk_sizes_pao, blk_sizes_pri
|
||||
LOGICAL :: found
|
||||
REAL(dp), DIMENSION(:), POINTER :: H_evals
|
||||
REAL(dp), DIMENSION(:, :), POINTER :: A, block_H0, block_N, block_N_inv, &
|
||||
block_X, H, H_evecs, V0
|
||||
TYPE(dbcsr_iterator_type) :: iter
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
|
||||
CALL dbcsr_get_info(pao%matrix_Y, row_blk_size=blk_sizes_pri, col_blk_size=blk_sizes_pao)
|
||||
|
||||
!$OMP PARALLEL DEFAULT(NONE) SHARED(pao,qs_env,blk_sizes_pri,blk_sizes_pao) &
|
||||
!$OMP PRIVATE(iter,arow,acol,iatom,n,m,i,found) &
|
||||
!$OMP PRIVATE(block_X,block_H0,block_N,block_N_inv,A,H,H_evecs,H_evals,V0)
|
||||
CALL dbcsr_iterator_start(iter, pao%matrix_X)
|
||||
DO WHILE (dbcsr_iterator_blocks_left(iter))
|
||||
CALL dbcsr_iterator_next_block(iter, arow, acol, block_X)
|
||||
iatom = arow; CPASSERT(arow == acol)
|
||||
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_H0, row=iatom, col=iatom, block=block_H0, found=found)
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_N_diag, row=iatom, col=iatom, block=block_N, found=found)
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_N_inv_diag, row=iatom, col=iatom, block=block_N_inv, found=found)
|
||||
CPASSERT(ASSOCIATED(block_H0) .AND. ASSOCIATED(block_N) .AND. ASSOCIATED(block_N_inv))
|
||||
|
||||
n = blk_sizes_pri(iatom) ! size of primary basis
|
||||
m = blk_sizes_pao(iatom) ! size of pao basis
|
||||
|
||||
ALLOCATE (V0(n, n))
|
||||
CALL pao_guess_initial_potential(qs_env, iatom, V0)
|
||||
|
||||
! construct H
|
||||
ALLOCATE (H(n, n))
|
||||
H = MATMUL(MATMUL(block_N, block_H0 + V0), block_N) ! transform into orthonormal basis
|
||||
|
||||
! diagonalize H
|
||||
ALLOCATE (H_evecs(n, n), H_evals(n))
|
||||
H_evecs = H
|
||||
CALL diamat_all(H_evecs, H_evals)
|
||||
|
||||
! use first m eigenvectors as initial guess
|
||||
ALLOCATE (A(n, m))
|
||||
A = MATMUL(block_N_inv, H_evecs(:, 1:m))
|
||||
|
||||
! normalize vectors
|
||||
DO i = 1, m
|
||||
A(:, i) = A(:, i)/NORM2(A(:, i))
|
||||
END DO
|
||||
|
||||
block_X = RESHAPE(A, (/n*m, 1/))
|
||||
DEALLOCATE (H, V0, A, H_evecs, H_evals)
|
||||
|
||||
END DO
|
||||
CALL dbcsr_iterator_stop(iter)
|
||||
!$OMP END PARALLEL
|
||||
|
||||
CALL timestop(handle)
|
||||
|
||||
END SUBROUTINE pao_param_initguess_equi
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Takes current matrix_X and calculates the matrices A and B.
|
||||
!> \param pao ...
|
||||
!> \param qs_env ...
|
||||
!> \param ls_scf_env ...
|
||||
!> \param gradient ...
|
||||
!> \param penalty ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_AB_equi(pao, qs_env, ls_scf_env, gradient, penalty)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
TYPE(ls_scf_env_type), TARGET :: ls_scf_env
|
||||
LOGICAL, INTENT(IN) :: gradient
|
||||
REAL(dp), INTENT(INOUT), OPTIONAL :: penalty
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_AB_equi'
|
||||
|
||||
INTEGER :: acol, arow, handle, i, iatom, j, k, m, n
|
||||
LOGICAL :: found
|
||||
REAL(dp) :: denom, w
|
||||
REAL(dp), DIMENSION(:), POINTER :: ANNA_evals
|
||||
REAL(dp), DIMENSION(:, :), POINTER :: ANNA, ANNA_evecs, ANNA_inv, block_A, &
|
||||
block_B, block_G, block_Ma, block_Mb, &
|
||||
block_N, block_X, D, G, M1, M2, M3, &
|
||||
M4, M5, NN
|
||||
TYPE(dbcsr_iterator_type) :: iter
|
||||
TYPE(dbcsr_type) :: matrix_Ma, matrix_Mb
|
||||
TYPE(ls_mstruct_type), POINTER :: ls_mstruct
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
ls_mstruct => ls_scf_env%ls_mstruct
|
||||
|
||||
IF (gradient) THEN
|
||||
CALL pao_calc_grad_lnv_wrt_AB(qs_env, ls_scf_env, matrix_Ma, matrix_Mb)
|
||||
END IF
|
||||
|
||||
!$OMP PARALLEL DEFAULT(NONE) SHARED(pao,ls_mstruct,matrix_Ma,matrix_Mb,gradient,penalty) &
|
||||
!$OMP PRIVATE(iter,arow,acol,iatom,found,n,m,w,i,j,k,denom) &
|
||||
!$OMP PRIVATE(NN,ANNA,ANNA_evals,ANNA_evecs,ANNA_inv,D,G,M1,M2,M3,M4,M5) &
|
||||
!$OMP PRIVATE(block_X,block_A,block_B,block_N,block_Ma, block_Mb, block_G)
|
||||
CALL dbcsr_iterator_start(iter, pao%matrix_X)
|
||||
DO WHILE (dbcsr_iterator_blocks_left(iter))
|
||||
CALL dbcsr_iterator_next_block(iter, arow, acol, block_X)
|
||||
iatom = arow; CPASSERT(arow == acol)
|
||||
CALL dbcsr_get_block_p(matrix=ls_mstruct%matrix_A, row=iatom, col=iatom, block=block_A, found=found)
|
||||
CPASSERT(ASSOCIATED(block_A))
|
||||
CALL dbcsr_get_block_p(matrix=ls_mstruct%matrix_B, row=iatom, col=iatom, block=block_B, found=found)
|
||||
CPASSERT(ASSOCIATED(block_B))
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_N_diag, row=iatom, col=iatom, block=block_N, found=found)
|
||||
CPASSERT(ASSOCIATED(block_N))
|
||||
|
||||
n = SIZE(block_A, 1) ! size of primary basis
|
||||
m = SIZE(block_A, 2) ! size of pao basis
|
||||
block_A = RESHAPE(block_X, (/n, m/))
|
||||
|
||||
! restrain pao basis vectors to unit norm
|
||||
IF (PRESENT(penalty)) THEN
|
||||
DO i = 1, m
|
||||
w = 1.0_dp - SUM(block_A(:, i)**2)
|
||||
penalty = penalty + pao%penalty_strength*w**2
|
||||
END DO
|
||||
END IF
|
||||
|
||||
ALLOCATE (NN(n, n), ANNA(m, m))
|
||||
NN = MATMUL(block_N, block_N) ! it's actually S^{-1}
|
||||
ANNA = MATMUL(MATMUL(TRANSPOSE(block_A), NN), block_A)
|
||||
|
||||
! diagonalize ANNA
|
||||
ALLOCATE (ANNA_evecs(m, m), ANNA_evals(m))
|
||||
ANNA_evecs(:, :) = ANNA
|
||||
CALL diamat_all(ANNA_evecs, ANNA_evals)
|
||||
IF (MINVAL(ABS(ANNA_evals)) < 1e-10_dp) CPABORT("PAO basis singualar.")
|
||||
|
||||
! build ANNA_inv
|
||||
ALLOCATE (ANNA_inv(m, m))
|
||||
ANNA_inv(:, :) = 0.0_dp
|
||||
DO k = 1, m
|
||||
w = 1.0_dp/ANNA_evals(k)
|
||||
DO i = 1, m
|
||||
DO j = 1, m
|
||||
ANNA_inv(i, j) = ANNA_inv(i, j) + w*ANNA_evecs(i, k)*ANNA_evecs(j, k)
|
||||
END DO
|
||||
END DO
|
||||
END DO
|
||||
|
||||
!B = 1/S * A * 1/(A^T 1/S A)
|
||||
block_B = MATMUL(MATMUL(NN, block_A), ANNA_inv)
|
||||
|
||||
! TURNING POINT (if calc grad) ------------------------------------------
|
||||
IF (gradient) THEN
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_G, row=iatom, col=iatom, block=block_G, found=found)
|
||||
CPASSERT(ASSOCIATED(block_G))
|
||||
CALL dbcsr_get_block_p(matrix=matrix_Ma, row=iatom, col=iatom, block=block_Ma, found=found)
|
||||
CALL dbcsr_get_block_p(matrix=matrix_Mb, row=iatom, col=iatom, block=block_Mb, found=found)
|
||||
! don't check ASSOCIATED(block_M), it might have been filtered out.
|
||||
|
||||
ALLOCATE (G(n, m))
|
||||
G(:, :) = 0.0_dp
|
||||
|
||||
IF (PRESENT(penalty)) THEN
|
||||
DO i = 1, m
|
||||
w = 1.0_dp - SUM(block_A(:, i)**2)
|
||||
G(:, i) = -4.0_dp*pao%penalty_strength*w*block_A(:, i)
|
||||
END DO
|
||||
END IF
|
||||
|
||||
IF (ASSOCIATED(block_Ma)) THEN
|
||||
G = G + block_Ma
|
||||
END IF
|
||||
|
||||
IF (ASSOCIATED(block_Mb)) THEN
|
||||
G = G + MATMUL(MATMUL(NN, block_Mb), ANNA_inv)
|
||||
|
||||
! calculate derivatives dAA_inv/ dAA
|
||||
ALLOCATE (D(m, m), M1(m, m), M2(m, m), M3(m, m), M4(m, m), M5(m, m))
|
||||
|
||||
DO i = 1, m
|
||||
DO j = 1, m
|
||||
denom = ANNA_evals(i) - ANNA_evals(j)
|
||||
IF (i == j) THEN
|
||||
D(i, i) = -1.0_dp/ANNA_evals(i)**2 ! diagonal elements
|
||||
ELSE IF (ABS(denom) > 1e-10_dp) THEN
|
||||
D(i, j) = (1.0_dp/ANNA_evals(i) - 1.0_dp/ANNA_evals(j))/denom
|
||||
ELSE
|
||||
D(i, j) = -1.0_dp ! limit according to L'Hospital's rule
|
||||
END IF
|
||||
END DO
|
||||
END DO
|
||||
|
||||
M1 = MATMUL(MATMUL(TRANSPOSE(block_A), NN), block_Mb)
|
||||
M2 = MATMUL(MATMUL(TRANSPOSE(ANNA_evecs), M1), ANNA_evecs)
|
||||
M3 = M2*D ! Hadamard product
|
||||
M4 = MATMUL(MATMUL(ANNA_evecs, M3), TRANSPOSE(ANNA_evecs))
|
||||
M5 = 0.5_dp*(M4 + TRANSPOSE(M4))
|
||||
G = G + 2.0_dp*MATMUL(MATMUL(NN, block_A), M5)
|
||||
|
||||
DEALLOCATE (D, M1, M2, M3, M4, M5)
|
||||
END IF
|
||||
|
||||
block_G = RESHAPE(G, (/n*m, 1/))
|
||||
DEALLOCATE (G)
|
||||
END IF
|
||||
|
||||
DEALLOCATE (NN, ANNA, ANNA_evecs, ANNA_evals, ANNA_inv)
|
||||
END DO
|
||||
CALL dbcsr_iterator_stop(iter)
|
||||
!$OMP END PARALLEL
|
||||
|
||||
IF (gradient) THEN
|
||||
CALL dbcsr_release(matrix_Ma)
|
||||
CALL dbcsr_release(matrix_Mb)
|
||||
END IF
|
||||
|
||||
CALL timestop(handle)
|
||||
|
||||
END SUBROUTINE pao_calc_AB_equi
|
||||
|
||||
END MODULE pao_param_equi
|
||||
|
|
@ -15,8 +15,11 @@ MODULE pao_param_exp
|
|||
dbcsr_create, dbcsr_get_block_p, dbcsr_get_info, dbcsr_iterator_blocks_left, &
|
||||
dbcsr_iterator_next_block, dbcsr_iterator_start, dbcsr_iterator_stop, dbcsr_iterator_type, &
|
||||
dbcsr_p_type, dbcsr_release, dbcsr_reserve_diag_blocks, dbcsr_set, dbcsr_type
|
||||
USE dm_ls_scf_types, ONLY: ls_scf_env_type
|
||||
USE kinds, ONLY: dp
|
||||
USE mathlib, ONLY: diamat_all
|
||||
USE pao_param_methods, ONLY: pao_calc_AB_from_U,&
|
||||
pao_calc_grad_lnv_wrt_U
|
||||
USE pao_potentials, ONLY: pao_guess_initial_potential
|
||||
USE pao_types, ONLY: pao_env_type
|
||||
USE qs_environment_types, ONLY: get_qs_env,&
|
||||
|
|
@ -31,7 +34,7 @@ MODULE pao_param_exp
|
|||
|
||||
CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'pao_param_exp'
|
||||
|
||||
PUBLIC :: pao_param_init_exp, pao_param_finalize_exp, pao_calc_U_exp
|
||||
PUBLIC :: pao_param_init_exp, pao_param_finalize_exp, pao_calc_AB_exp
|
||||
PUBLIC :: pao_param_count_exp, pao_param_initguess_exp
|
||||
|
||||
CONTAINS
|
||||
|
|
@ -156,14 +159,54 @@ CONTAINS
|
|||
|
||||
END SUBROUTINE pao_param_initguess_exp
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Takes current matrix_X and calculates the matrices A and B.
|
||||
!> \param pao ...
|
||||
!> \param qs_env ...
|
||||
!> \param ls_scf_env ...
|
||||
!> \param gradient ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_AB_exp(pao, qs_env, ls_scf_env, gradient)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
TYPE(ls_scf_env_type), TARGET :: ls_scf_env
|
||||
LOGICAL, INTENT(IN) :: gradient
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_AB_exp'
|
||||
|
||||
INTEGER :: handle
|
||||
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
|
||||
TYPE(dbcsr_type) :: matrix_M, matrix_U
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
CALL get_qs_env(qs_env, matrix_s=matrix_s)
|
||||
CALL dbcsr_create(matrix_U, matrix_type="N", dist=pao%diag_distribution, template=matrix_s(1)%matrix)
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_U)
|
||||
|
||||
!TODO: move this condition into pao_calc_U, use matrix_N as template
|
||||
IF (gradient) THEN
|
||||
CALL pao_calc_grad_lnv_wrt_U(qs_env, ls_scf_env, matrix_M)
|
||||
CALL pao_calc_U_exp(pao, matrix_U, matrix_M, pao%matrix_G)
|
||||
CALL dbcsr_release(matrix_M)
|
||||
ELSE
|
||||
CALL pao_calc_U_exp(pao, matrix_U)
|
||||
END IF
|
||||
|
||||
CALL pao_calc_AB_from_U(pao, qs_env, ls_scf_env, matrix_U)
|
||||
CALL dbcsr_release(matrix_U)
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_calc_AB_exp
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Calculate new matrix U and optionally its gradient G
|
||||
!> \param pao ...
|
||||
!> \param matrix_U ...
|
||||
!> \param matrix_M ...
|
||||
!> \param matrix_G ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_U_exp(pao, matrix_M, matrix_G)
|
||||
SUBROUTINE pao_calc_U_exp(pao, matrix_U, matrix_M, matrix_G)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(dbcsr_type) :: matrix_U
|
||||
TYPE(dbcsr_type), OPTIONAL :: matrix_M, matrix_G
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_U_exp'
|
||||
|
|
@ -184,7 +227,7 @@ CONTAINS
|
|||
|
||||
CALL dbcsr_get_info(pao%matrix_Y, row_blk_size=blk_sizes_pri, col_blk_size=blk_sizes_pao)
|
||||
|
||||
!$OMP PARALLEL DEFAULT(NONE) SHARED(pao,matrix_M,matrix_G,blk_sizes_pri,blk_sizes_pao) &
|
||||
!$OMP PARALLEL DEFAULT(NONE) SHARED(pao,matrix_U,matrix_M,matrix_G,blk_sizes_pri,blk_sizes_pao) &
|
||||
!$OMP PRIVATE(iter,arow,acol,iatom,N,M,nparams,i,j,k,found) &
|
||||
!$OMP PRIVATE(block_X,block_U,block_U0,block_X_full,evals,evecs) &
|
||||
!$OMP PRIVATE(block_M,block_G,block_D,block_tmp,block_G_full,denom)
|
||||
|
|
@ -192,7 +235,7 @@ CONTAINS
|
|||
DO WHILE (dbcsr_iterator_blocks_left(iter))
|
||||
CALL dbcsr_iterator_next_block(iter, arow, acol, block_X)
|
||||
iatom = arow; CPASSERT(arow == acol)
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_U, row=iatom, col=iatom, block=block_U, found=found)
|
||||
CALL dbcsr_get_block_p(matrix=matrix_U, row=iatom, col=iatom, block=block_U, found=found)
|
||||
CPASSERT(ASSOCIATED(block_U))
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_U0, row=iatom, col=iatom, block=block_U0, found=found)
|
||||
CPASSERT(ASSOCIATED(block_U0))
|
||||
|
|
|
|||
|
|
@ -29,6 +29,8 @@ MODULE pao_param_gth
|
|||
mp_sum
|
||||
USE orbital_pointers, ONLY: init_orbital_pointers
|
||||
USE pao_param_fock, ONLY: pao_calc_U_block_fock
|
||||
USE pao_param_methods, ONLY: pao_calc_AB_from_U,&
|
||||
pao_calc_grad_lnv_wrt_U
|
||||
USE pao_potentials, ONLY: pao_calc_gaussian
|
||||
USE pao_types, ONLY: pao_env_type
|
||||
USE particle_types, ONLY: particle_type
|
||||
|
|
@ -43,7 +45,7 @@ MODULE pao_param_gth
|
|||
|
||||
PRIVATE
|
||||
|
||||
PUBLIC :: pao_param_init_gth, pao_param_finalize_gth, pao_calc_U_gth
|
||||
PUBLIC :: pao_param_init_gth, pao_param_finalize_gth, pao_calc_AB_gth
|
||||
PUBLIC :: pao_param_count_gth, pao_param_initguess_gth
|
||||
|
||||
CONTAINS
|
||||
|
|
@ -239,17 +241,59 @@ CONTAINS
|
|||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_param_gth_preconditioner
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Takes current matrix_X and calculates the matrices A and B.
|
||||
!> \param pao ...
|
||||
!> \param qs_env ...
|
||||
!> \param ls_scf_env ...
|
||||
!> \param gradient ...
|
||||
!> \param penalty ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_AB_gth(pao, qs_env, ls_scf_env, gradient, penalty)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
TYPE(ls_scf_env_type), TARGET :: ls_scf_env
|
||||
LOGICAL, INTENT(IN) :: gradient
|
||||
REAL(dp), INTENT(INOUT), OPTIONAL :: penalty
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_AB_gth'
|
||||
|
||||
INTEGER :: handle
|
||||
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
|
||||
TYPE(dbcsr_type) :: matrix_M, matrix_U
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
CALL get_qs_env(qs_env, matrix_s=matrix_s)
|
||||
CALL dbcsr_create(matrix_U, matrix_type="N", dist=pao%diag_distribution, template=matrix_s(1)%matrix)
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_U)
|
||||
|
||||
!TODO: move this condition into pao_calc_U, use matrix_N as template
|
||||
IF (gradient) THEN
|
||||
CALL pao_calc_grad_lnv_wrt_U(qs_env, ls_scf_env, matrix_M)
|
||||
CALL pao_calc_U_gth(pao, matrix_U, matrix_M, pao%matrix_G, penalty)
|
||||
CALL dbcsr_release(matrix_M)
|
||||
ELSE
|
||||
CALL pao_calc_U_gth(pao, matrix_U, penalty=penalty)
|
||||
END IF
|
||||
|
||||
CALL pao_calc_AB_from_U(pao, qs_env, ls_scf_env, matrix_U)
|
||||
CALL dbcsr_release(matrix_U)
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_calc_AB_gth
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Calculate new matrix U and optinally its gradient G
|
||||
!> \param pao ...
|
||||
!> \param penalty ...
|
||||
!> \param matrix_U ...
|
||||
!> \param matrix_M1 ...
|
||||
!> \param matrix_G ...
|
||||
!> \param penalty ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_U_gth(pao, penalty, matrix_M1, matrix_G)
|
||||
SUBROUTINE pao_calc_U_gth(pao, matrix_U, matrix_M1, matrix_G, penalty)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
REAL(dp), INTENT(INOUT), OPTIONAL :: penalty
|
||||
TYPE(dbcsr_type) :: matrix_U
|
||||
TYPE(dbcsr_type), OPTIONAL :: matrix_M1, matrix_G
|
||||
REAL(dp), INTENT(INOUT), OPTIONAL :: penalty
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_U_gth'
|
||||
|
||||
|
|
@ -289,7 +333,7 @@ CONTAINS
|
|||
CALL mp_sum(world_X, group) ! sync world view across MPI ranks
|
||||
|
||||
! loop over atoms
|
||||
CALL dbcsr_iterator_start(iter, pao%matrix_U)
|
||||
CALL dbcsr_iterator_start(iter, matrix_U)
|
||||
DO WHILE (dbcsr_iterator_blocks_left(iter))
|
||||
CALL dbcsr_iterator_next_block(iter, arow, acol, block_U)
|
||||
iatom = arow; CPASSERT(arow == acol)
|
||||
|
|
|
|||
|
|
@ -17,7 +17,8 @@ MODULE pao_param_linpot
|
|||
USE dbcsr_api, ONLY: &
|
||||
dbcsr_create, dbcsr_get_block_p, dbcsr_get_info, dbcsr_iterator_blocks_left, &
|
||||
dbcsr_iterator_next_block, dbcsr_iterator_start, dbcsr_iterator_stop, dbcsr_iterator_type, &
|
||||
dbcsr_release, dbcsr_reserve_diag_blocks, dbcsr_type
|
||||
dbcsr_p_type, dbcsr_release, dbcsr_reserve_diag_blocks, dbcsr_type
|
||||
USE dm_ls_scf_types, ONLY: ls_scf_env_type
|
||||
USE kinds, ONLY: dp
|
||||
USE machine, ONLY: m_flush
|
||||
USE mathlib, ONLY: diamat_all
|
||||
|
|
@ -33,6 +34,8 @@ MODULE pao_param_linpot
|
|||
linpot_rotinv_calc_terms,&
|
||||
linpot_rotinv_count_terms
|
||||
USE pao_param_fock, ONLY: pao_calc_U_block_fock
|
||||
USE pao_param_methods, ONLY: pao_calc_AB_from_U,&
|
||||
pao_calc_grad_lnv_wrt_U
|
||||
USE pao_potentials, ONLY: pao_guess_initial_potential
|
||||
USE pao_types, ONLY: pao_env_type
|
||||
USE particle_types, ONLY: particle_type
|
||||
|
|
@ -46,7 +49,7 @@ MODULE pao_param_linpot
|
|||
|
||||
PRIVATE
|
||||
|
||||
PUBLIC :: pao_param_init_linpot, pao_param_finalize_linpot, pao_calc_U_linpot
|
||||
PUBLIC :: pao_param_init_linpot, pao_param_finalize_linpot, pao_calc_AB_linpot
|
||||
PUBLIC :: pao_param_count_linpot, pao_param_initguess_linpot
|
||||
|
||||
CONTAINS
|
||||
|
|
@ -342,20 +345,64 @@ CONTAINS
|
|||
|
||||
END SUBROUTINE pao_param_count_linpot
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Takes current matrix_X and calculates the matrices A and B.
|
||||
!> \param pao ...
|
||||
!> \param qs_env ...
|
||||
!> \param ls_scf_env ...
|
||||
!> \param gradient ...
|
||||
!> \param penalty ...
|
||||
!> \param forces ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_AB_linpot(pao, qs_env, ls_scf_env, gradient, penalty, forces)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
TYPE(ls_scf_env_type), TARGET :: ls_scf_env
|
||||
LOGICAL, INTENT(IN) :: gradient
|
||||
REAL(dp), INTENT(INOUT), OPTIONAL :: penalty
|
||||
REAL(dp), DIMENSION(:, :), INTENT(INOUT), OPTIONAL :: forces
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_AB_linpot'
|
||||
|
||||
INTEGER :: handle
|
||||
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
|
||||
TYPE(dbcsr_type) :: matrix_M, matrix_U
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
CALL get_qs_env(qs_env, matrix_s=matrix_s)
|
||||
CALL dbcsr_create(matrix_U, matrix_type="N", dist=pao%diag_distribution, template=matrix_s(1)%matrix)
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_U)
|
||||
|
||||
!TODO: move this condition into pao_calc_U, use matrix_N as template
|
||||
IF (gradient) THEN
|
||||
CALL pao_calc_grad_lnv_wrt_U(qs_env, ls_scf_env, matrix_M)
|
||||
CALL pao_calc_U_linpot(pao, qs_env, matrix_U, matrix_M, pao%matrix_G, penalty, forces)
|
||||
CALL dbcsr_release(matrix_M)
|
||||
ELSE
|
||||
CALL pao_calc_U_linpot(pao, qs_env, matrix_U, penalty=penalty)
|
||||
END IF
|
||||
|
||||
CALL pao_calc_AB_from_U(pao, qs_env, ls_scf_env, matrix_U)
|
||||
CALL dbcsr_release(matrix_U)
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_calc_AB_linpot
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Calculate new matrix U and optinally its gradient G
|
||||
!> \param pao ...
|
||||
!> \param qs_env ...
|
||||
!> \param penalty ...
|
||||
!> \param matrix_U ...
|
||||
!> \param matrix_M ...
|
||||
!> \param matrix_G ...
|
||||
!> \param penalty ...
|
||||
!> \param forces ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_U_linpot(pao, qs_env, penalty, matrix_M, matrix_G, forces)
|
||||
SUBROUTINE pao_calc_U_linpot(pao, qs_env, matrix_U, matrix_M, matrix_G, penalty, forces)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
REAL(dp), INTENT(INOUT), OPTIONAL :: penalty
|
||||
TYPE(dbcsr_type) :: matrix_U
|
||||
TYPE(dbcsr_type), OPTIONAL :: matrix_M, matrix_G
|
||||
REAL(dp), INTENT(INOUT), OPTIONAL :: penalty
|
||||
REAL(dp), DIMENSION(:, :), INTENT(INOUT), OPTIONAL :: forces
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_U_linpot'
|
||||
|
|
@ -383,7 +430,7 @@ CONTAINS
|
|||
evals(:, :) = 0.0_dp
|
||||
gaps(:) = HUGE(1.0_dp)
|
||||
regu_energy = 0.0_dp
|
||||
CALL dbcsr_get_info(pao%matrix_U, group=group_handle)
|
||||
CALL dbcsr_get_info(matrix_U, group=group_handle)
|
||||
CALL group%set_handle(group_handle)
|
||||
|
||||
CALL dbcsr_iterator_start(iter, pao%matrix_X)
|
||||
|
|
@ -391,7 +438,7 @@ CONTAINS
|
|||
CALL dbcsr_iterator_next_block(iter, arow, acol, block_X)
|
||||
iatom = arow; CPASSERT(arow == acol)
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_R, row=iatom, col=iatom, block=block_R, found=found)
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_U, row=iatom, col=iatom, block=block_U, found=found)
|
||||
CALL dbcsr_get_block_p(matrix=matrix_U, row=iatom, col=iatom, block=block_U, found=found)
|
||||
CPASSERT(ASSOCIATED(block_R) .AND. ASSOCIATED(block_U))
|
||||
n = SIZE(block_U, 1)
|
||||
|
||||
|
|
@ -415,17 +462,18 @@ CONTAINS
|
|||
IF (PRESENT(penalty) .AND. nterms > 0) &
|
||||
regu_energy = regu_energy + DOT_PRODUCT(block_X(:, 1), MATMUL(block_R, block_X(:, 1)))
|
||||
|
||||
IF (.NOT. PRESENT(matrix_G) .AND. .NOT. PRESENT(matrix_G)) THEN
|
||||
CALL pao_calc_U_block_fock(pao, iatom=iatom, penalty=penalty, V=block_V, U=block_U, &
|
||||
gap=gaps(iatom), evals=evals(:, iatom))
|
||||
CALL pao_calc_U_block_fock(pao, iatom=iatom, penalty=penalty, V=block_V, U=block_U, &
|
||||
gap=gaps(iatom), evals=evals(:, iatom))
|
||||
|
||||
ELSE ! TURNING POINT (if calc grad) -------------------------------------------------------
|
||||
IF (PRESENT(matrix_G)) THEN ! TURNING POINT (if calc grad) --------------------------------
|
||||
CPASSERT(PRESENT(matrix_M))
|
||||
CALL dbcsr_get_block_p(matrix=matrix_M, row=iatom, col=iatom, block=block_M1, found=found)
|
||||
|
||||
! corner-cases: block_M1 might have been filtered out or there might be zero pao parameters
|
||||
IF (ASSOCIATED(block_M1) .AND. SIZE(block_V_terms) > 0) THEN
|
||||
ALLOCATE (vec_M2(n*n))
|
||||
block_M2(1:n, 1:n) => vec_M2(:) ! map vector into matrix
|
||||
!TODO: this 2nd call does double work. However, *sometimes* this branch is not taken.
|
||||
CALL pao_calc_U_block_fock(pao, iatom=iatom, penalty=penalty, V=block_V, U=block_U, &
|
||||
M1=block_M1, G=block_M2, gap=gaps(iatom), evals=evals(:, iatom))
|
||||
IF (MAXVAL(ABS(block_M2 - TRANSPOSE(block_M2))) > 1e-14_dp) &
|
||||
|
|
|
|||
444
src/pao_param_methods.F
Normal file
444
src/pao_param_methods.F
Normal file
|
|
@ -0,0 +1,444 @@
|
|||
!--------------------------------------------------------------------------------------------------!
|
||||
! CP2K: A general program to perform molecular dynamics simulations !
|
||||
! Copyright 2000-2023 CP2K developers group <https://cp2k.org> !
|
||||
! !
|
||||
! SPDX-License-Identifier: GPL-2.0-or-later !
|
||||
!--------------------------------------------------------------------------------------------------!
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Common routines for PAO parametrizations.
|
||||
!> \author Ole Schuett
|
||||
! **************************************************************************************************
|
||||
MODULE pao_param_methods
|
||||
USE cp_control_types, ONLY: dft_control_type
|
||||
USE cp_log_handling, ONLY: cp_to_string
|
||||
USE dbcsr_api, ONLY: &
|
||||
dbcsr_add, dbcsr_complete_redistribute, dbcsr_create, dbcsr_get_block_p, dbcsr_get_info, &
|
||||
dbcsr_iterator_blocks_left, dbcsr_iterator_next_block, dbcsr_iterator_start, &
|
||||
dbcsr_iterator_stop, dbcsr_iterator_type, dbcsr_multiply, dbcsr_p_type, dbcsr_release, &
|
||||
dbcsr_reserve_diag_blocks, dbcsr_scale, dbcsr_type
|
||||
USE dm_ls_scf_qs, ONLY: matrix_decluster
|
||||
USE dm_ls_scf_types, ONLY: ls_mstruct_type,&
|
||||
ls_scf_env_type
|
||||
USE kinds, ONLY: dp
|
||||
USE message_passing, ONLY: mp_comm_type,&
|
||||
mp_max
|
||||
USE pao_types, ONLY: pao_env_type
|
||||
USE qs_environment_types, ONLY: get_qs_env,&
|
||||
qs_environment_type
|
||||
USE qs_rho_types, ONLY: qs_rho_get,&
|
||||
qs_rho_type
|
||||
#include "./base/base_uses.f90"
|
||||
|
||||
IMPLICIT NONE
|
||||
|
||||
PRIVATE
|
||||
|
||||
CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'pao_param_methods'
|
||||
|
||||
PUBLIC :: pao_calc_grad_lnv_wrt_U, pao_calc_AB_from_U, pao_calc_grad_lnv_wrt_AB
|
||||
|
||||
CONTAINS
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Helper routine, calculates partial derivative dE/dU
|
||||
!> \param qs_env ...
|
||||
!> \param ls_scf_env ...
|
||||
!> \param matrix_M_diag the derivate wrt U, matrix uses pao%diag_distribution
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_grad_lnv_wrt_U(qs_env, ls_scf_env, matrix_M_diag)
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
TYPE(ls_scf_env_type), TARGET :: ls_scf_env
|
||||
TYPE(dbcsr_type) :: matrix_M_diag
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_grad_lnv_wrt_U', &
|
||||
routineP = moduleN//':'//routineN
|
||||
|
||||
INTEGER :: handle
|
||||
REAL(KIND=dp) :: filter_eps
|
||||
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
|
||||
TYPE(dbcsr_type) :: matrix_M, matrix_Ma, matrix_Mb, matrix_NM
|
||||
TYPE(ls_mstruct_type), POINTER :: ls_mstruct
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
|
||||
ls_mstruct => ls_scf_env%ls_mstruct
|
||||
pao => ls_scf_env%pao_env
|
||||
filter_eps = ls_scf_env%eps_filter
|
||||
CALL get_qs_env(qs_env, matrix_s=matrix_s)
|
||||
|
||||
CALL pao_calc_grad_lnv_wrt_AB(qs_env, ls_scf_env, matrix_Ma, matrix_Mb)
|
||||
|
||||
! Calculation uses distr. of matrix_s, afterwards we redistribute to pao%diag_distribution.
|
||||
CALL dbcsr_create(matrix_M, template=matrix_s(1)%matrix, matrix_type="N")
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_M)
|
||||
|
||||
CALL dbcsr_create(matrix_NM, template=ls_mstruct%matrix_A, matrix_type="N")
|
||||
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, pao%matrix_N_inv, matrix_Ma, &
|
||||
1.0_dp, matrix_NM, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, pao%matrix_N, matrix_Mb, &
|
||||
1.0_dp, matrix_NM, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "T", 1.0_dp, matrix_NM, pao%matrix_Y, &
|
||||
1.0_dp, matrix_M, filter_eps=filter_eps)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! redistribute using pao%diag_distribution
|
||||
CALL dbcsr_create(matrix_M_diag, &
|
||||
name="PAO matrix_M", &
|
||||
matrix_type="N", &
|
||||
dist=pao%diag_distribution, &
|
||||
template=matrix_s(1)%matrix)
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_M_diag)
|
||||
CALL dbcsr_complete_redistribute(matrix_M, matrix_M_diag)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! cleanup:
|
||||
CALL dbcsr_release(matrix_M)
|
||||
CALL dbcsr_release(matrix_Ma)
|
||||
CALL dbcsr_release(matrix_Mb)
|
||||
CALL dbcsr_release(matrix_NM)
|
||||
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_calc_grad_lnv_wrt_U
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Takes current matrix_X and calculates the matrices A and B.
|
||||
!> \param pao ...
|
||||
!> \param qs_env ...
|
||||
!> \param ls_scf_env ...
|
||||
!> \param matrix_U_diag ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_AB_from_U(pao, qs_env, ls_scf_env, matrix_U_diag)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
TYPE(ls_scf_env_type), TARGET :: ls_scf_env
|
||||
TYPE(dbcsr_type) :: matrix_U_diag
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_AB_from_U'
|
||||
|
||||
INTEGER :: acol, arow, handle, iatom
|
||||
LOGICAL :: found
|
||||
REAL(dp), DIMENSION(:, :), POINTER :: block_A, block_B, block_N, block_N_inv, &
|
||||
block_U, block_Y
|
||||
TYPE(dbcsr_iterator_type) :: iter
|
||||
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
|
||||
TYPE(dbcsr_type) :: matrix_U
|
||||
TYPE(ls_mstruct_type), POINTER :: ls_mstruct
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
CALL get_qs_env(qs_env, matrix_s=matrix_s)
|
||||
ls_mstruct => ls_scf_env%ls_mstruct
|
||||
|
||||
! --------------------------------------------------------------------------------------------
|
||||
! sanity check matrix U
|
||||
CALL pao_assert_unitary(pao, matrix_U_diag)
|
||||
|
||||
! --------------------------------------------------------------------------------------------
|
||||
! redistribute matrix_U_diag from diag_distribution to distribution of matrix_s
|
||||
CALL get_qs_env(qs_env, matrix_s=matrix_s)
|
||||
CALL dbcsr_create(matrix_U, matrix_type="N", template=matrix_s(1)%matrix)
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_U)
|
||||
CALL dbcsr_complete_redistribute(matrix_U_diag, matrix_U)
|
||||
|
||||
! --------------------------------------------------------------------------------------------
|
||||
! calculate matrix A and B from matrix U
|
||||
! Multiplying diagonal matrices is a local operation.
|
||||
! To take advantage of this we're using an iterator instead of calling dbcsr_multiply().
|
||||
!$OMP PARALLEL DEFAULT(NONE) SHARED(pao,ls_mstruct,matrix_U) &
|
||||
!$OMP PRIVATE(iter,arow,acol,iatom,block_U,block_Y,block_A,block_B,block_N,block_N_inv,found)
|
||||
CALL dbcsr_iterator_start(iter, matrix_U)
|
||||
DO WHILE (dbcsr_iterator_blocks_left(iter))
|
||||
CALL dbcsr_iterator_next_block(iter, arow, acol, block_U)
|
||||
iatom = arow; CPASSERT(arow == acol)
|
||||
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_Y, row=iatom, col=iatom, block=block_Y, found=found)
|
||||
CPASSERT(ASSOCIATED(block_Y))
|
||||
|
||||
CALL dbcsr_get_block_p(matrix=ls_mstruct%matrix_A, row=iatom, col=iatom, block=block_A, found=found)
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_N_inv, row=iatom, col=iatom, block=block_N_inv, found=found)
|
||||
CPASSERT(ASSOCIATED(block_A) .AND. ASSOCIATED(block_N_inv))
|
||||
|
||||
CALL dbcsr_get_block_p(matrix=ls_mstruct%matrix_B, row=iatom, col=iatom, block=block_B, found=found)
|
||||
CALL dbcsr_get_block_p(matrix=pao%matrix_N, row=iatom, col=iatom, block=block_N, found=found)
|
||||
CPASSERT(ASSOCIATED(block_B) .AND. ASSOCIATED(block_N))
|
||||
|
||||
block_A = MATMUL(MATMUL(block_N_inv, block_U), block_Y)
|
||||
block_B = MATMUL(MATMUL(block_N, block_U), block_Y)
|
||||
|
||||
END DO
|
||||
CALL dbcsr_iterator_stop(iter)
|
||||
!$OMP END PARALLEL
|
||||
|
||||
CALL dbcsr_release(matrix_U)
|
||||
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_calc_AB_from_U
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Debugging routine, check unitaryness of U
|
||||
!> \param pao ...
|
||||
!> \param matrix_U ...
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_assert_unitary(pao, matrix_U)
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(dbcsr_type) :: matrix_U
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_assert_unitary', &
|
||||
routineP = moduleN//':'//routineN
|
||||
|
||||
INTEGER :: acol, arow, group_handle, handle, i, &
|
||||
iatom, M, N
|
||||
INTEGER, DIMENSION(:), POINTER :: blk_sizes_pao, blk_sizes_pri
|
||||
REAL(dp) :: delta_max
|
||||
REAL(dp), DIMENSION(:, :), POINTER :: block_test, tmp1, tmp2
|
||||
TYPE(dbcsr_iterator_type) :: iter
|
||||
TYPE(mp_comm_type) :: group
|
||||
|
||||
IF (pao%check_unitary_tol < 0.0_dp) RETURN ! no checking
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
delta_max = 0.0_dp
|
||||
|
||||
CALL dbcsr_get_info(pao%matrix_Y, row_blk_size=blk_sizes_pri, col_blk_size=blk_sizes_pao)
|
||||
|
||||
!$OMP PARALLEL DEFAULT(NONE) SHARED(pao,matrix_U,blk_sizes_pri,blk_sizes_pao,delta_max) &
|
||||
!$OMP PRIVATE(iter,arow,acol,iatom,N,M,block_test,tmp1,tmp2)
|
||||
CALL dbcsr_iterator_start(iter, matrix_U)
|
||||
DO WHILE (dbcsr_iterator_blocks_left(iter))
|
||||
CALL dbcsr_iterator_next_block(iter, arow, acol, block_test)
|
||||
iatom = arow; CPASSERT(arow == acol)
|
||||
N = blk_sizes_pri(iatom) ! size of primary basis
|
||||
M = blk_sizes_pao(iatom) ! size of pao basis
|
||||
|
||||
! we only need the upper left "PAO-corner" to be unitary
|
||||
ALLOCATE (tmp1(N, M), tmp2(M, M))
|
||||
tmp1 = block_test(:, 1:M)
|
||||
tmp2 = MATMUL(TRANSPOSE(tmp1), tmp1)
|
||||
DO i = 1, M
|
||||
tmp2(i, i) = tmp2(i, i) - 1.0_dp
|
||||
END DO
|
||||
|
||||
!$OMP ATOMIC
|
||||
delta_max = MAX(delta_max, MAXVAL(ABS(tmp2)))
|
||||
|
||||
DEALLOCATE (tmp1, tmp2)
|
||||
END DO
|
||||
CALL dbcsr_iterator_stop(iter)
|
||||
!$OMP END PARALLEL
|
||||
|
||||
CALL dbcsr_get_info(matrix_U, group=group_handle)
|
||||
CALL group%set_handle(group_handle)
|
||||
|
||||
CALL mp_max(delta_max, group)
|
||||
IF (pao%iw > 0) WRITE (pao%iw, *) 'PAO| checked unitaryness, max delta:', delta_max
|
||||
IF (delta_max > pao%check_unitary_tol) &
|
||||
CPABORT("Found bad unitaryness:"//cp_to_string(delta_max))
|
||||
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_assert_unitary
|
||||
|
||||
! **************************************************************************************************
|
||||
!> \brief Helper routine, calculates partial derivative dE/dA and dE/dB.
|
||||
!> As energy functional serves the definition by LNV (Li, Nunes, Vanderbilt).
|
||||
!> \param qs_env ...
|
||||
!> \param ls_scf_env ...
|
||||
!> \param matrix_Ma the derivate wrt A, matrix uses s_matrix-distribution.
|
||||
!> \param matrix_Mb the derivate wrt B, matrix uses s_matrix-distribution.
|
||||
! **************************************************************************************************
|
||||
SUBROUTINE pao_calc_grad_lnv_wrt_AB(qs_env, ls_scf_env, matrix_Ma, matrix_Mb)
|
||||
TYPE(qs_environment_type), POINTER :: qs_env
|
||||
TYPE(ls_scf_env_type), TARGET :: ls_scf_env
|
||||
TYPE(dbcsr_type) :: matrix_Ma, matrix_Mb
|
||||
|
||||
CHARACTER(len=*), PARAMETER :: routineN = 'pao_calc_grad_lnv_wrt_AB', &
|
||||
routineP = moduleN//':'//routineN
|
||||
|
||||
INTEGER :: handle, nspin
|
||||
INTEGER, DIMENSION(:), POINTER :: pao_blk_sizes
|
||||
REAL(KIND=dp) :: filter_eps
|
||||
TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s, rho_ao
|
||||
TYPE(dbcsr_type) :: matrix_HB, matrix_HPS, matrix_M, matrix_M1, matrix_M1_dc, matrix_M2, &
|
||||
matrix_M2_dc, matrix_M3, matrix_M3_dc, matrix_PA, matrix_PH, matrix_PHP, matrix_PSP, &
|
||||
matrix_SB, matrix_SP
|
||||
TYPE(dft_control_type), POINTER :: dft_control
|
||||
TYPE(ls_mstruct_type), POINTER :: ls_mstruct
|
||||
TYPE(pao_env_type), POINTER :: pao
|
||||
TYPE(qs_rho_type), POINTER :: rho
|
||||
|
||||
CALL timeset(routineN, handle)
|
||||
|
||||
ls_mstruct => ls_scf_env%ls_mstruct
|
||||
pao => ls_scf_env%pao_env
|
||||
|
||||
CALL get_qs_env(qs_env, &
|
||||
rho=rho, &
|
||||
matrix_ks=matrix_ks, &
|
||||
matrix_s=matrix_s, &
|
||||
dft_control=dft_control)
|
||||
CALL qs_rho_get(rho, rho_ao=rho_ao)
|
||||
nspin = dft_control%nspins
|
||||
filter_eps = ls_scf_env%eps_filter
|
||||
|
||||
CALL dbcsr_get_info(ls_mstruct%matrix_A, col_blk_size=pao_blk_sizes)
|
||||
|
||||
IF (nspin /= 1) CPABORT("open shell not yet implemented")
|
||||
!TODO: handle openshell case properly
|
||||
|
||||
! Notation according to equation (4.6) on page 50 from:
|
||||
! https://dx.doi.org/10.3929%2Fethz-a-010819495
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! calculate need products in pao basis
|
||||
CALL dbcsr_create(matrix_PH, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, ls_scf_env%matrix_p(1), ls_scf_env%matrix_ks(1), &
|
||||
0.0_dp, matrix_PH, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_PHP, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_PH, ls_scf_env%matrix_p(1), &
|
||||
0.0_dp, matrix_PHP, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_SP, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, ls_scf_env%matrix_s, ls_scf_env%matrix_p(1), &
|
||||
0.0_dp, matrix_SP, filter_eps=filter_eps)
|
||||
|
||||
IF (nspin == 1) CALL dbcsr_scale(matrix_SP, 0.5_dp)
|
||||
|
||||
CALL dbcsr_create(matrix_HPS, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "T", 1.0_dp, ls_scf_env%matrix_ks(1), matrix_SP, &
|
||||
0.0_dp, matrix_HPS, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_create(matrix_PSP, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, ls_scf_env%matrix_p(1), matrix_SP, &
|
||||
0.0_dp, matrix_PSP, filter_eps=filter_eps)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! M1 = dE_lnv / dP_pao
|
||||
CALL dbcsr_create(matrix_M1, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
|
||||
CALL dbcsr_multiply("N", "T", 3.0_dp, ls_scf_env%matrix_ks(1), matrix_SP, &
|
||||
1.0_dp, matrix_M1, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "N", 3.0_dp, matrix_SP, ls_scf_env%matrix_ks(1), &
|
||||
1.0_dp, matrix_M1, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "T", -2.0_dp, matrix_HPS, matrix_SP, &
|
||||
1.0_dp, matrix_M1, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "N", -2.0_dp, matrix_SP, matrix_HPS, &
|
||||
1.0_dp, matrix_M1, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "T", -2.0_dp, matrix_SP, matrix_HPS, &
|
||||
1.0_dp, matrix_M1, filter_eps=filter_eps)
|
||||
|
||||
! reverse possible molecular clustering
|
||||
CALL dbcsr_create(matrix_M1_dc, &
|
||||
template=matrix_s(1)%matrix, &
|
||||
row_blk_size=pao_blk_sizes, &
|
||||
col_blk_size=pao_blk_sizes)
|
||||
CALL matrix_decluster(matrix_M1_dc, matrix_M1, ls_mstruct)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! M2 = dE_lnv / dH
|
||||
CALL dbcsr_create(matrix_M2, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
|
||||
CALL dbcsr_add(matrix_M2, matrix_PSP, 1.0_dp, 3.0_dp)
|
||||
|
||||
CALL dbcsr_multiply("N", "N", -2.0_dp, matrix_PSP, matrix_SP, &
|
||||
1.0_dp, matrix_M2, filter_eps=filter_eps)
|
||||
|
||||
! reverse possible molecular clustering
|
||||
CALL dbcsr_create(matrix_M2_dc, &
|
||||
template=matrix_s(1)%matrix, &
|
||||
row_blk_size=pao_blk_sizes, &
|
||||
col_blk_size=pao_blk_sizes)
|
||||
CALL matrix_decluster(matrix_M2_dc, matrix_M2, ls_mstruct)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! M3 = dE_lnv / dS
|
||||
CALL dbcsr_create(matrix_M3, template=ls_scf_env%matrix_s, matrix_type="N")
|
||||
|
||||
CALL dbcsr_add(matrix_M3, matrix_PHP, 1.0_dp, 3.0_dp)
|
||||
|
||||
CALL dbcsr_multiply("N", "N", -2.0_dp, matrix_PHP, matrix_SP, &
|
||||
1.0_dp, matrix_M3, filter_eps=filter_eps)
|
||||
|
||||
CALL dbcsr_multiply("N", "T", -2.0_dp, matrix_PSP, matrix_PH, &
|
||||
1.0_dp, matrix_M3, filter_eps=filter_eps)
|
||||
|
||||
! reverse possible molecular clustering
|
||||
CALL dbcsr_create(matrix_M3_dc, &
|
||||
template=matrix_s(1)%matrix, &
|
||||
row_blk_size=pao_blk_sizes, &
|
||||
col_blk_size=pao_blk_sizes)
|
||||
CALL matrix_decluster(matrix_M3_dc, matrix_M3, ls_mstruct)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! assemble Ma and Mb
|
||||
! matrix_Ma = dE_lnv / dA = P * A * M1
|
||||
! matrix_Mb = dE_lnv / dB = H * B * M2 + S * B * M3
|
||||
CALL dbcsr_create(matrix_Ma, template=ls_mstruct%matrix_A, matrix_type="N")
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_Ma)
|
||||
CALL dbcsr_create(matrix_Mb, template=ls_mstruct%matrix_B, matrix_type="N")
|
||||
CALL dbcsr_reserve_diag_blocks(matrix_Mb)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! combine M1 with matrices from primary basis
|
||||
CALL dbcsr_create(matrix_PA, template=ls_mstruct%matrix_A, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, rho_ao(1)%matrix, ls_mstruct%matrix_A, &
|
||||
0.0_dp, matrix_PA, filter_eps=filter_eps)
|
||||
|
||||
! matrix_Ma = P * A * M1
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_PA, matrix_M1_dc, &
|
||||
0.0_dp, matrix_Ma, filter_eps=filter_eps)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! combine M2 with matrices from primary basis
|
||||
CALL dbcsr_create(matrix_HB, template=ls_mstruct%matrix_B, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_ks(1)%matrix, ls_mstruct%matrix_B, &
|
||||
0.0_dp, matrix_HB, filter_eps=filter_eps)
|
||||
|
||||
! matrix_Mb = H * B * M2
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_HB, matrix_M2_dc, &
|
||||
0.0_dp, matrix_Mb, filter_eps=filter_eps)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! combine M3 with matrices from primary basis
|
||||
CALL dbcsr_create(matrix_SB, template=ls_mstruct%matrix_B, matrix_type="N")
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_s(1)%matrix, ls_mstruct%matrix_B, &
|
||||
0.0_dp, matrix_SB, filter_eps=filter_eps)
|
||||
|
||||
IF (nspin == 1) CALL dbcsr_scale(matrix_SB, 0.5_dp)
|
||||
|
||||
! matrix_Mb += S * B * M3
|
||||
CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_SB, matrix_M3_dc, &
|
||||
1.0_dp, matrix_Mb, filter_eps=filter_eps)
|
||||
|
||||
IF (nspin == 1) CALL dbcsr_scale(matrix_Ma, 2.0_dp)
|
||||
IF (nspin == 1) CALL dbcsr_scale(matrix_Mb, 2.0_dp)
|
||||
|
||||
!---------------------------------------------------------------------------
|
||||
! cleanup: TODO release matrices as early as possible
|
||||
CALL dbcsr_release(matrix_PH)
|
||||
CALL dbcsr_release(matrix_PHP)
|
||||
CALL dbcsr_release(matrix_SP)
|
||||
CALL dbcsr_release(matrix_HPS)
|
||||
CALL dbcsr_release(matrix_PSP)
|
||||
CALL dbcsr_release(matrix_M)
|
||||
CALL dbcsr_release(matrix_M1)
|
||||
CALL dbcsr_release(matrix_M2)
|
||||
CALL dbcsr_release(matrix_M3)
|
||||
CALL dbcsr_release(matrix_M1_dc)
|
||||
CALL dbcsr_release(matrix_M2_dc)
|
||||
CALL dbcsr_release(matrix_M3_dc)
|
||||
CALL dbcsr_release(matrix_PA)
|
||||
CALL dbcsr_release(matrix_HB)
|
||||
CALL dbcsr_release(matrix_SB)
|
||||
|
||||
CALL timestop(handle)
|
||||
END SUBROUTINE pao_calc_grad_lnv_wrt_AB
|
||||
|
||||
END MODULE pao_param_methods
|
||||
|
|
@ -100,14 +100,14 @@ MODULE pao_types
|
|||
!> \var constants_ready set when stuff, which does not depend of atomic positions is ready
|
||||
!> \var need_initial_scf set when the initial density matrix is not self-consistend
|
||||
!> \var matrix_X parameters of pao basis, which eventually determine matrix_U. Uses diag_distribution.
|
||||
!> \var matrix_U roation matrix derived from matrix_X. Uses diag_distribution.
|
||||
!> \var matrix_U0 constant pre-rotation which serves as initial guess for exp-parametrization. Uses diag_distribution.
|
||||
!> \var matrix_H0 Diagonal blocks of core hamiltonian, uses diag_distribution
|
||||
!> \var matrix_Y selector matrix which translates between primary and pao basis.
|
||||
!> basically a block diagonal "rectangular identity matrix". Uses s_matrix-distribution.
|
||||
!> \var matrix_N diagonal matrix filled with 1/sqrt(S) from primary overlap matrix. Uses s_matrix-distribution.
|
||||
!> \var matrix_N_inv diagonal matrix filled with sqrt(S) from primary overlap matrix. Uses s_matrix-distribution.
|
||||
!> \var matrix_N_diag copy of matrix_N using diag_distribution
|
||||
!> \var matrix_N_inv diagonal matrix filled with sqrt(S) from primary overlap matrix. Uses s_matrix-distribution.
|
||||
!> \var matrix_N_inv_diag copy of matrix_N_inv using diag_distribution
|
||||
!> \var matrix_X_orig copy made of matrix_X at beginning of optimization cylce, used for mixing. Uses diag_distribution.
|
||||
!> \var matrix_G derivative of pao-energy wrt to matrix_X, ie. the pao-gradient. Uses diag_distribution.
|
||||
!> \var matrix_G_prev copy of gradient from previous step, used for conjugate gradient method. Uses diag_distribution.
|
||||
|
|
@ -180,13 +180,13 @@ MODULE pao_types
|
|||
|
||||
! matrices
|
||||
TYPE(dbcsr_type) :: matrix_X
|
||||
TYPE(dbcsr_type) :: matrix_U
|
||||
TYPE(dbcsr_type) :: matrix_U0
|
||||
TYPE(dbcsr_type) :: matrix_H0
|
||||
TYPE(dbcsr_type) :: matrix_Y
|
||||
TYPE(dbcsr_type) :: matrix_N
|
||||
TYPE(dbcsr_type) :: matrix_N_inv
|
||||
TYPE(dbcsr_type) :: matrix_N_diag
|
||||
TYPE(dbcsr_type) :: matrix_N_inv
|
||||
TYPE(dbcsr_type) :: matrix_N_inv_diag
|
||||
TYPE(dbcsr_type) :: matrix_X_orig
|
||||
TYPE(dbcsr_type) :: matrix_G
|
||||
TYPE(dbcsr_type) :: matrix_G_prev
|
||||
|
|
@ -222,8 +222,9 @@ CONTAINS
|
|||
CALL dbcsr_release(pao%matrix_X)
|
||||
CALL dbcsr_release(pao%matrix_Y)
|
||||
CALL dbcsr_release(pao%matrix_N)
|
||||
CALL dbcsr_release(pao%matrix_N_inv)
|
||||
CALL dbcsr_release(pao%matrix_N_diag)
|
||||
CALL dbcsr_release(pao%matrix_N_inv)
|
||||
CALL dbcsr_release(pao%matrix_N_inv_diag)
|
||||
CALL dbcsr_release(pao%matrix_H0)
|
||||
|
||||
DEALLOCATE (pao%ml_training_set)
|
||||
|
|
|
|||
58
tests/QS/regtest-pao-1/H2_pao_equi.inp
Normal file
58
tests/QS/regtest-pao-1/H2_pao_equi.inp
Normal file
|
|
@ -0,0 +1,58 @@
|
|||
&FORCE_EVAL
|
||||
METHOD Quickstep
|
||||
&DFT
|
||||
BASIS_SET_FILE_NAME BASIS_MOLOPT
|
||||
POTENTIAL_FILE_NAME GTH_POTENTIALS
|
||||
&QS
|
||||
LS_SCF
|
||||
&END QS
|
||||
&POISSON
|
||||
PERIODIC NONE
|
||||
PSOLVER MT
|
||||
&END POISSON
|
||||
&LS_SCF
|
||||
MAX_SCF 25
|
||||
EPS_SCF 1.0E-8
|
||||
S_PRECONDITIONER NONE
|
||||
REPORT_ALL_SPARSITIES OFF
|
||||
PURIFICATION_METHOD TRS4
|
||||
EXTRAPOLATION_ORDER 1
|
||||
&PAO
|
||||
MAX_PAO 500
|
||||
EPS_PAO 1.0E-5
|
||||
PARAMETERIZATION EQUIVARIANT
|
||||
CHECK_UNITARY_TOL 1.0E-10
|
||||
&LINE_SEARCH
|
||||
METHOD GOLD
|
||||
&END LINE_SEARCH
|
||||
&END PAO
|
||||
&END
|
||||
&XC
|
||||
&XC_FUNCTIONAL PBE
|
||||
&END XC_FUNCTIONAL
|
||||
&END XC
|
||||
&END DFT
|
||||
&SUBSYS
|
||||
&CELL
|
||||
ABC 4.0 4.0 4.0
|
||||
PERIODIC NONE
|
||||
&END CELL
|
||||
&COORD
|
||||
H 0.72 0.0 0.0
|
||||
H 0.0 0.0 0.0
|
||||
&END COORD
|
||||
&KIND H
|
||||
BASIS_SET DZVP-MOLOPT-SR-GTH
|
||||
POTENTIAL GTH-PBE
|
||||
PAO_BASIS_SIZE 1
|
||||
&END KIND
|
||||
&TOPOLOGY
|
||||
&CENTER_COORDINATES
|
||||
&END
|
||||
&END
|
||||
&END SUBSYS
|
||||
&END FORCE_EVAL
|
||||
&GLOBAL
|
||||
PROJECT H2_pao_equi
|
||||
RUN_TYPE ENERGY
|
||||
&END GLOBAL
|
||||
59
tests/QS/regtest-pao-1/H2_pao_equi_checkgrad.inp
Normal file
59
tests/QS/regtest-pao-1/H2_pao_equi_checkgrad.inp
Normal file
|
|
@ -0,0 +1,59 @@
|
|||
&FORCE_EVAL
|
||||
METHOD Quickstep
|
||||
&DFT
|
||||
BASIS_SET_FILE_NAME BASIS_MOLOPT
|
||||
POTENTIAL_FILE_NAME GTH_POTENTIALS
|
||||
&QS
|
||||
LS_SCF
|
||||
&END QS
|
||||
&POISSON
|
||||
PERIODIC NONE
|
||||
PSOLVER MT
|
||||
&END POISSON
|
||||
&LS_SCF
|
||||
MAX_SCF 25
|
||||
EPS_SCF 1.0E-8
|
||||
S_PRECONDITIONER NONE
|
||||
REPORT_ALL_SPARSITIES OFF
|
||||
PURIFICATION_METHOD TRS4
|
||||
EXTRAPOLATION_ORDER 1
|
||||
&PAO
|
||||
MAX_PAO 1
|
||||
EPS_PAO 1.0E-5
|
||||
PARAMETERIZATION EQUIVARIANT
|
||||
NUM_GRADIENT_ORDER 6
|
||||
CHECK_GRADIENT_TOL 1.0E-6
|
||||
&LINE_SEARCH
|
||||
METHOD 3PNT
|
||||
&END LINE_SEARCH
|
||||
&END PAO
|
||||
&END
|
||||
&XC
|
||||
&XC_FUNCTIONAL PBE
|
||||
&END XC_FUNCTIONAL
|
||||
&END XC
|
||||
&END DFT
|
||||
&SUBSYS
|
||||
&CELL
|
||||
ABC 4.0 4.0 4.0
|
||||
PERIODIC NONE
|
||||
&END CELL
|
||||
&COORD
|
||||
H 0.72 0.0 0.0
|
||||
H 0.0 0.0 0.0
|
||||
&END COORD
|
||||
&KIND H
|
||||
BASIS_SET DZVP-MOLOPT-SR-GTH
|
||||
POTENTIAL GTH-PBE
|
||||
PAO_BASIS_SIZE 1
|
||||
&END KIND
|
||||
&TOPOLOGY
|
||||
&CENTER_COORDINATES
|
||||
&END
|
||||
&END
|
||||
&END SUBSYS
|
||||
&END FORCE_EVAL
|
||||
&GLOBAL
|
||||
PROJECT H2_pao_equi_checkgrad
|
||||
RUN_TYPE ENERGY
|
||||
&END GLOBAL
|
||||
|
|
@ -4,6 +4,7 @@ H2_pao_exp.inp 11 1e-07
|
|||
H2_pao_fock.inp 11 1e-07 -1.160920966911983
|
||||
H2_pao_rotinv.inp 11 1e-07 -1.160920966911983
|
||||
H2_pao_gth.inp 11 1e-07 -1.160920966911983
|
||||
H2_pao_equi.inp 11 1e-07 -1.160920966911983
|
||||
#
|
||||
H2_pao_MD.inp 2 1e-07 -0.116085621247E+01
|
||||
#
|
||||
|
|
@ -11,6 +12,7 @@ H2_pao_exp_checkgrad.inp 0
|
|||
H2_pao_fock_checkgrad.inp 0
|
||||
H2_pao_rotinv_checkgrad.inp 0
|
||||
H2_pao_gth_checkgrad.inp 0
|
||||
H2_pao_equi_checkgrad.inp 0
|
||||
#
|
||||
H2_pao_fock_checkforces.inp 0
|
||||
H2_pao_rotinv_checkforces.inp 0
|
||||
|
|
|
|||
71
tests/QS/regtest-pao-2/H2O_pao_equi.inp
Normal file
71
tests/QS/regtest-pao-2/H2O_pao_equi.inp
Normal file
|
|
@ -0,0 +1,71 @@
|
|||
&FORCE_EVAL
|
||||
METHOD Quickstep
|
||||
&DFT
|
||||
BASIS_SET_FILE_NAME BASIS_MOLOPT
|
||||
POTENTIAL_FILE_NAME GTH_POTENTIALS
|
||||
&QS
|
||||
LS_SCF
|
||||
&END QS
|
||||
&POISSON
|
||||
PERIODIC NONE
|
||||
PSOLVER MT
|
||||
&END POISSON
|
||||
&LS_SCF
|
||||
MAX_SCF 25
|
||||
EPS_SCF 1.0E-8
|
||||
EPS_FILTER 1.0E-8
|
||||
S_PRECONDITIONER NONE
|
||||
REPORT_ALL_SPARSITIES OFF
|
||||
PURIFICATION_METHOD TRS4
|
||||
EXTRAPOLATION_ORDER 1
|
||||
&PAO
|
||||
PREOPT_DM_FILE H2O_ref_LS_DM_SPIN_1_RESTART.dm
|
||||
EPS_PAO 1.0E-6
|
||||
MAX_PAO 3000
|
||||
PARAMETERIZATION EQUIVARIANT
|
||||
&LINE_SEARCH
|
||||
METHOD ADAPT
|
||||
&END LINE_SEARCH
|
||||
&PRINT
|
||||
&RESTART ON
|
||||
&END RESTART
|
||||
&ATOM_INFO ON
|
||||
&END ATOM_INFO
|
||||
&END PRINT
|
||||
&END PAO
|
||||
&END
|
||||
&XC
|
||||
&XC_FUNCTIONAL PBE
|
||||
&END XC_FUNCTIONAL
|
||||
&END XC
|
||||
&END DFT
|
||||
&SUBSYS
|
||||
&CELL
|
||||
ABC 4.0 4.0 4.0
|
||||
PERIODIC NONE
|
||||
&END CELL
|
||||
&COORD
|
||||
O 2.6116774290 4.1629472392 4.1629480502
|
||||
H 1.8450021896 4.5850871565 4.5850871241
|
||||
H 2.2373131663 3.5277354616 3.5277353471
|
||||
&END COORD
|
||||
&KIND H
|
||||
BASIS_SET DZVP-MOLOPT-SR-GTH
|
||||
POTENTIAL GTH-PBE
|
||||
PAO_BASIS_SIZE 4
|
||||
&END KIND
|
||||
&KIND O
|
||||
BASIS_SET DZVP-MOLOPT-SR-GTH
|
||||
POTENTIAL GTH-PBE
|
||||
PAO_BASIS_SIZE 4
|
||||
&END KIND
|
||||
&TOPOLOGY
|
||||
&CENTER_COORDINATES
|
||||
&END
|
||||
&END
|
||||
&END SUBSYS
|
||||
&END FORCE_EVAL
|
||||
&GLOBAL
|
||||
PROJECT H2O_pao_equi
|
||||
RUN_TYPE ENERGY
|
||||
&END GLOBAL
|
||||
69
tests/QS/regtest-pao-2/H2O_pao_equi_hybrid.inp
Normal file
69
tests/QS/regtest-pao-2/H2O_pao_equi_hybrid.inp
Normal file
|
|
@ -0,0 +1,69 @@
|
|||
&FORCE_EVAL
|
||||
METHOD Quickstep
|
||||
&DFT
|
||||
BASIS_SET_FILE_NAME BASIS_MOLOPT
|
||||
POTENTIAL_FILE_NAME GTH_POTENTIALS
|
||||
&QS
|
||||
LS_SCF
|
||||
&END QS
|
||||
&POISSON
|
||||
PERIODIC NONE
|
||||
PSOLVER MT
|
||||
&END POISSON
|
||||
&LS_SCF
|
||||
MAX_SCF 25
|
||||
EPS_SCF 1.0E-8
|
||||
EPS_FILTER 1.0E-8
|
||||
S_PRECONDITIONER NONE
|
||||
REPORT_ALL_SPARSITIES OFF
|
||||
PURIFICATION_METHOD TRS4
|
||||
EXTRAPOLATION_ORDER 1
|
||||
&PAO
|
||||
PREOPT_DM_FILE H2O_ref_LS_DM_SPIN_1_RESTART.dm
|
||||
EPS_PAO 1.0E-6
|
||||
MAX_PAO 3000
|
||||
PARAMETERIZATION EQUIVARIANT
|
||||
&LINE_SEARCH
|
||||
METHOD ADAPT
|
||||
&END LINE_SEARCH
|
||||
&PRINT
|
||||
&RESTART ON
|
||||
&END RESTART
|
||||
&END PRINT
|
||||
&END PAO
|
||||
&END
|
||||
&XC
|
||||
&XC_FUNCTIONAL PBE
|
||||
&END XC_FUNCTIONAL
|
||||
&END XC
|
||||
&END DFT
|
||||
&SUBSYS
|
||||
&CELL
|
||||
ABC 4.0 4.0 4.0
|
||||
PERIODIC NONE
|
||||
&END CELL
|
||||
&COORD
|
||||
O 2.6116774290 4.1629472392 4.1629480502
|
||||
H 1.8450021896 4.5850871565 4.5850871241
|
||||
H 2.2373131663 3.5277354616 3.5277353471
|
||||
&END COORD
|
||||
&KIND H
|
||||
BASIS_SET DZVP-MOLOPT-SR-GTH
|
||||
POTENTIAL GTH-PBE
|
||||
PAO_BASIS_SIZE 4
|
||||
&END KIND
|
||||
&KIND O
|
||||
BASIS_SET DZVP-MOLOPT-SR-GTH
|
||||
POTENTIAL GTH-PBE
|
||||
! disable pao for Oxygen
|
||||
&END KIND
|
||||
&TOPOLOGY
|
||||
&CENTER_COORDINATES
|
||||
&END
|
||||
&END
|
||||
&END SUBSYS
|
||||
&END FORCE_EVAL
|
||||
&GLOBAL
|
||||
PROJECT H2O_pao_equi_hybrid
|
||||
RUN_TYPE ENERGY
|
||||
&END GLOBAL
|
||||
|
|
@ -26,6 +26,10 @@
|
|||
&LINE_SEARCH
|
||||
METHOD ADAPT
|
||||
&END LINE_SEARCH
|
||||
&PRINT
|
||||
&ATOM_INFO
|
||||
&END ATOM_INFO
|
||||
&END PRINT
|
||||
&END PAO
|
||||
&END
|
||||
&XC
|
||||
|
|
@ -60,6 +64,6 @@
|
|||
&END SUBSYS
|
||||
&END FORCE_EVAL
|
||||
&GLOBAL
|
||||
PROJECT H2O_pao_fock
|
||||
PROJECT H2O_pao_minimal
|
||||
RUN_TYPE ENERGY
|
||||
&END GLOBAL
|
||||
|
|
|
|||
|
|
@ -5,6 +5,7 @@ H2O_pao_fock.inp 11 1e-08 -
|
|||
H2O_pao_rotinv.inp 11 1e-06 -17.202700812326320
|
||||
H2O_pao_eq_prim.inp 11 3e-09 -17.202700812326320
|
||||
H2O_pao_gth.inp 11 9e-05 -17.202700812326320
|
||||
H2O_pao_equi.inp 11 1e-08 -17.202700812326320
|
||||
#
|
||||
# molecular blocking and single precision
|
||||
#
|
||||
|
|
@ -19,6 +20,7 @@ H2O_pao_exp_hybrid.inp 11 1e-08 -
|
|||
H2O_pao_fock_hybrid.inp 11 1e-08 -17.202700812326320
|
||||
H2O_pao_rotinv_hybrid.inp 11 1e-06 -17.202700812326320
|
||||
H2O_pao_gth_hybrid.inp 11 9e-05 -17.202700812326320
|
||||
H2O_pao_equi_hybrid.inp 11 1e-08 -17.202700812326320
|
||||
#
|
||||
# rotated system, energy changes slightly due to grid and cell
|
||||
#
|
||||
|
|
|
|||
53
tests/QS/regtest-pao-3/He_pao_initguess_equi.inp
Normal file
53
tests/QS/regtest-pao-3/He_pao_initguess_equi.inp
Normal file
|
|
@ -0,0 +1,53 @@
|
|||
&FORCE_EVAL
|
||||
METHOD Quickstep
|
||||
&DFT
|
||||
BASIS_SET_FILE_NAME BASIS_MOLOPT
|
||||
POTENTIAL_FILE_NAME GTH_POTENTIALS
|
||||
&QS
|
||||
LS_SCF
|
||||
&END QS
|
||||
&POISSON
|
||||
PERIODIC NONE
|
||||
PSOLVER MT
|
||||
&END POISSON
|
||||
&LS_SCF
|
||||
MAX_SCF 25
|
||||
EPS_SCF 1.0E-8
|
||||
EPS_FILTER 1.0E-8
|
||||
S_PRECONDITIONER NONE
|
||||
REPORT_ALL_SPARSITIES OFF
|
||||
PURIFICATION_METHOD TRS4
|
||||
EXTRAPOLATION_ORDER 1
|
||||
&PAO
|
||||
MAX_PAO 0
|
||||
PARAMETERIZATION EQUIVARIANT
|
||||
&END PAO
|
||||
&END
|
||||
&XC
|
||||
&XC_FUNCTIONAL PBE
|
||||
&END XC_FUNCTIONAL
|
||||
&END XC
|
||||
&END DFT
|
||||
&SUBSYS
|
||||
&CELL
|
||||
ABC 4.0 4.0 4.0
|
||||
PERIODIC NONE
|
||||
&END CELL
|
||||
&COORD
|
||||
He 2.0 2.0 2.0
|
||||
&END COORD
|
||||
&KIND He
|
||||
BASIS_SET DZVP-MOLOPT-SR-GTH
|
||||
POTENTIAL GTH-PBE
|
||||
PAO_BASIS_SIZE 1
|
||||
&END KIND
|
||||
&TOPOLOGY
|
||||
&CENTER_COORDINATES
|
||||
&END
|
||||
&END
|
||||
&END SUBSYS
|
||||
&END FORCE_EVAL
|
||||
&GLOBAL
|
||||
PROJECT He_pao_initguess_equi
|
||||
RUN_TYPE ENERGY
|
||||
&END GLOBAL
|
||||
|
|
@ -4,4 +4,5 @@ He_pao_initguess_exp.inp 11 5e-4
|
|||
He_pao_initguess_fock.inp 11 5e-4 -2.889334869153459
|
||||
He_pao_initguess_rotinv.inp 11 5e-4 -2.889334869153459
|
||||
He_pao_initguess_gth.inp 11 5e-4 -2.889334869153459
|
||||
He_pao_initguess_equi.inp 11 5e-4 -2.889334869153459
|
||||
#EOF
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue