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

xTB pairpotential (A. Hehn) (#1136)

Browse source 
* xTB pairpotential (A. Hehn)

* Adjust regtest

* Adjust regtests

* adjust test
This commit is contained in:
Juerg Hutter 2020-10-21 14:51:22 +02:00 • committed by GitHub
parent ca80fb4e41
commit fe04ab922d
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
16 changed files with 812 additions and 31 deletions
Download patch file Download diff file Expand all files Collapse all files

9
src/cp_control_types.F
View file

@ -16,6 +16,8 @@ MODULE cp_control_types
USE kinds, ONLY: default_path_length,&
default_string_length,&
dp
USE pair_potential_types, ONLY: pair_potential_p_release,&
pair_potential_p_type
USE qs_cdft_types, ONLY: cdft_control_create,&
cdft_control_release,&
cdft_control_type
@ -113,6 +115,7 @@ MODULE cp_control_types
REAL(KIND=dp) :: kxr, kx2
!
LOGICAL :: xb_interaction
LOGICAL :: do_nonbonded
LOGICAL :: coulomb_interaction
LOGICAL :: tb3_interaction
LOGICAL :: check_atomic_charges
@ -124,6 +127,8 @@ MODULE cp_control_types
INTEGER, DIMENSION(:, :), POINTER :: kab_types
INTEGER :: kab_nval
REAL, DIMENSION(:), POINTER :: kab_vals
!
TYPE(pair_potential_p_type), POINTER :: nonbonded
END TYPE xtb_control_type
! **************************************************************************************************
@ -1052,6 +1057,7 @@ CONTAINS
NULLIFY (xtb_control%kab_param)
NULLIFY (xtb_control%kab_vals)
NULLIFY (xtb_control%kab_types)
NULLIFY (xtb_control%nonbonded)
END SUBROUTINE xtb_control_create
@ -1072,6 +1078,9 @@ CONTAINS
IF (ASSOCIATED(xtb_control%kab_types)) THEN
DEALLOCATE (xtb_control%kab_types)
END IF
IF (ASSOCIATED(xtb_control%nonbonded)) THEN
CALL pair_potential_p_release(xtb_control%nonbonded)
END IF
DEALLOCATE (xtb_control)
END IF
END SUBROUTINE xtb_control_release

26
src/cp_control_utils.F
View file

@ -26,6 +26,7 @@ MODULE cp_control_utils
cp_print_key_unit_nr
USE cp_units, ONLY: cp_unit_from_cp2k,&
cp_unit_to_cp2k
USE force_fields_input, ONLY: read_gp_section
USE input_constants, ONLY: &
constant_env, custom_env, do_admm_basis_projection, do_admm_blocked_projection, &
do_admm_blocking_purify_full, do_admm_charge_constrained_projection, &
@ -54,6 +55,7 @@ MODULE cp_control_utils
USE kinds, ONLY: default_path_length,&
default_string_length,&
dp
USE pair_potential_types, ONLY: pair_potential_reallocate
USE periodic_table, ONLY: get_ptable_info
USE qs_cdft_utils, ONLY: read_cdft_control_section
USE string_utilities, ONLY: uppercase
@ -613,21 +615,21 @@ CONTAINS
CHARACTER(LEN=default_string_length), &
DIMENSION(:), POINTER :: clist
INTEGER :: handle, itmp, j, jj, k, n_rep, n_var, &
ngauss, nrep
ngauss, ngp, nrep
INTEGER, DIMENSION(:), POINTER :: tmplist
LOGICAL :: explicit, was_present
REAL(dp) :: tmp, tmpsqrt, value
REAL(dp), POINTER :: scal(:)
TYPE(section_vals_type), POINTER :: cdft_control_section, ddapc_restraint_section, &
dftb_parameter, dftb_section, lri_optbas_section, mull_section, s2_restraint_section, &
se_section, xtb_parameter, xtb_section
dftb_parameter, dftb_section, genpot_section, lri_optbas_section, mull_section, &
nonbonded_section, s2_restraint_section, se_section, xtb_parameter, xtb_section
CALL timeset(routineN, handle)
was_present = .FALSE.
NULLIFY (mull_section, ddapc_restraint_section, s2_restraint_section, &
se_section, dftb_section, xtb_section, dftb_parameter, xtb_parameter, lri_optbas_section, &
cdft_control_section)
cdft_control_section, genpot_section)
mull_section => section_vals_get_subs_vals(qs_section, "MULLIKEN_RESTRAINT")
ddapc_restraint_section => section_vals_get_subs_vals(qs_section, "DDAPC_RESTRAINT")
@ -639,6 +641,8 @@ CONTAINS
xtb_parameter => section_vals_get_subs_vals(xtb_section, "PARAMETER")
lri_optbas_section => section_vals_get_subs_vals(qs_section, "OPTIMIZE_LRI_BASIS")
cdft_control_section => section_vals_get_subs_vals(qs_section, "CDFT")
nonbonded_section => section_vals_get_subs_vals(xtb_section, "NONBONDED")
genpot_section => section_vals_get_subs_vals(nonbonded_section, "GENPOT")
! Setup all defaults values and overwrite input parameters
! EPS_DEFAULT should set the target accuracy in the total energy (~per electron) or a closely related value
@ -1106,6 +1110,18 @@ CONTAINS
CALL section_vals_val_get(xtb_parameter, "HALOGEN_BINDING", r_vals=scal)
qs_control%xtb_control%kxr = scal(1)
qs_control%xtb_control%kx2 = scal(2)
! NONBONDED interactions
CALL section_vals_val_get(xtb_section, "DO_NONBONDED", &
l_val=qs_control%xtb_control%do_nonbonded)
CALL section_vals_get(nonbonded_section, explicit=explicit)
IF (explicit .AND. qs_control%xtb_control%do_nonbonded) THEN
CALL section_vals_get(genpot_section, explicit=explicit, n_repetition=ngp)
IF (explicit) THEN
CALL pair_potential_reallocate(qs_control%xtb_control%nonbonded, 1, ngp, gp=.TRUE.)
CALL read_gp_section(qs_control%xtb_control%nonbonded, genpot_section, 0)
END IF
END IF !nonbonded
! XB_radius
CALL section_vals_val_get(xtb_parameter, "XB_RADIUS", r_val=qs_control%xtb_control%xb_radius)
! Kab
@ -1557,6 +1573,8 @@ CONTAINS
"xTB| Halogen interaction potential", xtb_control%xb_interaction
WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
"xTB| Halogen interaction potential cutoff radius", xtb_control%xb_radius
WRITE (UNIT=output_unit, FMT="(T2,A,T71,L10)") &
"xTB| Nonbonded interactions", xtb_control%do_nonbonded
WRITE (UNIT=output_unit, FMT="(T2,A,T31,A50)") &
"xTB| D3 Dispersion: Parameter file", ADJUSTR(TRIM(xtb_control%dispersion_parameter_file))
WRITE (UNIT=output_unit, FMT="(T2,A,T51,3F10.3)") &

9
src/input_cp2k_dft.F
View file

@ -1140,6 +1140,15 @@ CONTAINS
CALL section_add_keyword(print_key, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, __LOCATION__, &
name="sab_xtb_nonbond", &
description="Activates the printing of nonbonded interaction from (xTB) "// &
"neighbor lists ", &
default_l_val=.FALSE., &
lone_keyword_l_val=.TRUE.)
CALL section_add_keyword(print_key, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, __LOCATION__, &
name="soo_list", &
description="Activates the printing of RI orbital-orbital "// &

46
src/input_cp2k_tb.F
View file

@ -12,7 +12,6 @@
!> \author fawzi
! **************************************************************************************************
MODULE input_cp2k_tb
!qxtb citation
USE bibliography, ONLY: Elstner1998,&
Grimme2017,&
Hu2007,&
@ -22,6 +21,7 @@ MODULE input_cp2k_tb
USE input_constants, ONLY: dispersion_d3,&
dispersion_uff,&
slater
USE input_cp2k_mm, ONLY: create_GENPOT_section
USE input_keyword_types, ONLY: keyword_create,&
keyword_release,&
keyword_type
@ -142,6 +142,10 @@ CONTAINS
CALL section_add_subsection(section, subsection)
CALL section_release(subsection)
CALL create_xtb_nonbonded_section(subsection)
CALL section_add_subsection(section, subsection)
CALL section_release(subsection)
NULLIFY (keyword)
CALL keyword_create(keyword, __LOCATION__, name="DO_EWALD", &
description="Use Ewald type method instead of direct sum for Coulomb interaction", &
@ -167,6 +171,13 @@ CONTAINS
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, __LOCATION__, name="DO_NONBONDED", &
description="Controls the computation of real-space "// &
"(short-range) nonbonded interactions as correction to xTB.", &
usage="DO_NONBONDED T", default_l_val=.FALSE., lone_keyword_l_val=.TRUE.)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, __LOCATION__, name="COULOMB_INTERACTION", &
description="Use Coulomb interaction terms (electrostatics + TB3); for debug only", &
usage="COULOMB_INTERACTION T", default_l_val=.TRUE., lone_keyword_l_val=.TRUE.)
@ -395,7 +406,40 @@ CONTAINS
CALL keyword_release(keyword)
END SUBROUTINE create_xtb_parameter_section
! **************************************************************************************************
!> \brief ...
!> \param section ...
! **************************************************************************************************
SUBROUTINE create_xtb_nonbonded_section(section)
TYPE(section_type), POINTER :: section
TYPE(keyword_type), POINTER :: keyword
TYPE(section_type), POINTER :: subsection
CPASSERT(.NOT. ASSOCIATED(section))
CALL section_create(section, __LOCATION__, name="NONBONDED", &
description="This section specifies the input parameters for NON-BONDED interactions.", &
n_keywords=1, n_subsections=0, repeats=.FALSE.)
NULLIFY (subsection)
CALL create_GENPOT_section(subsection)
CALL section_add_subsection(section, subsection)
CALL section_release(subsection)
CALL keyword_create(keyword, __LOCATION__, name="DX", &
description="Parameter used for computing the derivative with the Ridders' method.", &
usage="DX <REAL>", default_r_val=0.1_dp, unit_str="bohr")
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, __LOCATION__, name="ERROR_LIMIT", &
description="Checks that the error in computing the derivative is not larger than "// &
"the value set; in case error is larger a warning message is printed.", &
usage="ERROR_LIMIT <REAL>", default_r_val=1.0E-12_dp)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
END SUBROUTINE create_xtb_nonbonded_section
! **************************************************************************************************
!> \brief Creates the &ATOM_APRAMETER section
!> \param section the section to create

3
src/pair_potential.F
View file

@ -62,7 +62,8 @@ MODULE pair_potential
INTEGER, PARAMETER, PRIVATE :: MAX_POINTS = 2000000
PUBLIC :: spline_nonbond_control, &
get_nonbond_storage
get_nonbond_storage, &
init_genpot
CONTAINS

4
src/qs_energy_types.F
View file

@ -70,7 +70,9 @@ MODULE qs_energy_types
kinetic, & !total kinetic energy [rk]
s_square, &
surf_dipole, &
embed_corr ! correction for embedding potential
embed_corr, & ! correction for embedding potential
xtb_xb_inter, & ! correction for halogen bonding within GFN1-xTB
xtb_nonbonded ! correction for nonbonded interactions within GFN1-xTB
REAL(KIND=dp), DIMENSION(:), POINTER :: ddapc_restraint
END TYPE qs_energy_type

6
src/qs_environment_types.F
View file

@ -364,6 +364,7 @@ CONTAINS
!> \param sab_tbe ...
!> \param sab_core ...
!> \param sab_xb ...
!> \param sab_xtb_nonbond ...
!> \param sab_almo ...
!> \param sab_kp ...
!> \param particle_set ...
@ -506,7 +507,7 @@ CONTAINS
SUBROUTINE get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, &
dft_control, mos, mos_aux_fit, sab_orb, sab_aux_fit, sab_aux_fit_asymm, sab_aux_fit_vs_orb, &
sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, &
sab_se, sab_tbe, sab_core, sab_xb, sab_almo, sab_kp, particle_set, energy, force, &
sab_se, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_almo, sab_kp, particle_set, energy, force, &
matrix_h, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, &
matrix_h_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_RI_aux_kp, &
matrix_ks_aux_fit, matrix_ks_aux_fit_im, matrix_ks_aux_fit_dft, matrix_ks_aux_fit_hfx, &
@ -547,7 +548,7 @@ CONTAINS
LOGICAL, OPTIONAL :: qmmm, qmmm_periodic
TYPE(neighbor_list_set_p_type), DIMENSION(:), OPTIONAL, POINTER :: sac_ae, sac_ppl, sac_lri, &
sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_tbe, sab_core, sab_xb, &
sab_almo, sab_kp
sab_xtb_nonbond, sab_almo, sab_kp
TYPE(particle_type), DIMENSION(:), OPTIONAL, &
POINTER :: particle_set
TYPE(qs_energy_type), OPTIONAL, POINTER :: energy
@ -830,6 +831,7 @@ CONTAINS
sab_tbe=sab_tbe, &
sab_core=sab_core, &
sab_xb=sab_xb, &
sab_xtb_nonbond=sab_xtb_nonbond, &
sab_almo=sab_almo, &
sab_kp=sab_kp, &
task_list=task_list, &

17
src/qs_ks_types.F
View file

@ -116,6 +116,7 @@ MODULE qs_ks_types
!> Ewald terms (DFTB)
!> \param sab_core neighbor lists for the calculation of the core interactions
!> \param sab_xb neighbor lists for the calculation of the XB interactions in xTB
!> \param sab_xtb_nonbond neighbor lists for the calculation of the nonbonded interactions in xTB
!> \param sab_all neighbor lists for the calculation of the matrix element of
!> non-symmetric operators
!> \param sab_vdw: neighbor lists for the calculation of dispersion interactions
@ -194,6 +195,7 @@ MODULE qs_ks_types
sab_tbe => Null(), &
sab_core => Null(), &
sab_xb => Null(), &
sab_xtb_nonbond => Null(), &
sab_all => Null(), &
sab_aux_fit => Null(), &
sab_aux_fit_vs_orb => Null(), &
@ -304,6 +306,7 @@ CONTAINS
!> \param sab_tbe ...
!> \param sab_core ...
!> \param sab_xb ...
!> \param sab_xtb_nonbond ...
!> \param sab_vdw ...
!> \param sab_scp ...
!> \param sab_almo ...
@ -357,7 +360,7 @@ CONTAINS
neighbor_list_id, &
sab_orb, sab_aux_fit, sab_aux_fit_asymm, &
sab_aux_fit_vs_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, &
sab_se, sab_tbe, sab_core, sab_xb, sab_vdw, sab_scp, sab_almo, sab_kp, &
sab_se, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, &
task_list, task_list_aux_fit, task_list_soft, &
kpoints, do_kpoints, &
atomic_kind_set, qs_kind_set, cell, cell_ref, use_ref_cell, &
@ -385,7 +388,8 @@ CONTAINS
INTEGER, OPTIONAL :: neighbor_list_id
TYPE(neighbor_list_set_p_type), DIMENSION(:), OPTIONAL, POINTER :: sab_orb, sab_aux_fit, &
sab_aux_fit_asymm, sab_aux_fit_vs_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, &
sap_oce, sab_lrc, sab_se, sab_tbe, sab_core, sab_xb, sab_vdw, sab_scp, sab_almo, sab_kp
sap_oce, sab_lrc, sab_se, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_vdw, sab_scp, &
sab_almo, sab_kp
TYPE(task_list_type), OPTIONAL, POINTER :: task_list, task_list_aux_fit, &
task_list_soft
TYPE(kpoint_type), OPTIONAL, POINTER :: kpoints
@ -486,6 +490,7 @@ CONTAINS
IF (PRESENT(sab_tbe)) sab_tbe => ks_env%sab_tbe
IF (PRESENT(sab_core)) sab_core => ks_env%sab_core
IF (PRESENT(sab_xb)) sab_xb => ks_env%sab_xb
IF (PRESENT(sab_xtb_nonbond)) sab_xtb_nonbond => ks_env%sab_xtb_nonbond
IF (PRESENT(sab_almo)) sab_almo => ks_env%sab_almo
IF (PRESENT(sab_kp)) sab_kp => ks_env%sab_kp
IF (PRESENT(dft_control)) dft_control => ks_env%dft_control
@ -583,6 +588,7 @@ CONTAINS
!> \param sab_tbe ...
!> \param sab_core ...
!> \param sab_xb ...
!> \param sab_xtb_nonbond ...
!> \param sab_vdw ...
!> \param sab_scp ...
!> \param sab_almo ...
@ -612,7 +618,7 @@ CONTAINS
kpoints, &
sab_orb, sab_aux_fit, sab_aux_fit_asymm, &
sab_aux_fit_vs_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, &
sab_se, sab_tbe, sab_core, sab_xb, sab_vdw, sab_scp, sab_almo, sab_kp, &
sab_se, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, &
task_list, task_list_aux_fit, task_list_soft, &
subsys, dft_control, dbcsr_dist, distribution_2d, pw_env, &
para_env, blacs_env, admm_dm)
@ -635,7 +641,8 @@ CONTAINS
TYPE(kpoint_type), OPTIONAL, POINTER :: kpoints
TYPE(neighbor_list_set_p_type), DIMENSION(:), OPTIONAL, POINTER :: sab_orb, sab_aux_fit, &
sab_aux_fit_asymm, sab_aux_fit_vs_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, &
sap_oce, sab_lrc, sab_se, sab_tbe, sab_core, sab_xb, sab_vdw, sab_scp, sab_almo, sab_kp
sap_oce, sab_lrc, sab_se, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_vdw, sab_scp, &
sab_almo, sab_kp
TYPE(task_list_type), OPTIONAL, POINTER :: task_list, task_list_aux_fit, &
task_list_soft
TYPE(qs_subsys_type), OPTIONAL, POINTER :: subsys
@ -713,6 +720,7 @@ CONTAINS
IF (PRESENT(sab_tbe)) ks_env%sab_tbe => sab_tbe
IF (PRESENT(sab_core)) ks_env%sab_core => sab_core
IF (PRESENT(sab_xb)) ks_env%sab_xb => sab_xb
IF (PRESENT(sab_xtb_nonbond)) ks_env%sab_xtb_nonbond => sab_xtb_nonbond
IF (PRESENT(sab_almo)) ks_env%sab_almo => sab_almo
IF (PRESENT(sab_kp)) ks_env%sab_kp => sab_kp
IF (PRESENT(admm_dm)) ks_env%admm_dm => admm_dm
@ -860,6 +868,7 @@ CONTAINS
CALL release_neighbor_list_sets(ks_env%sab_tbe)
CALL release_neighbor_list_sets(ks_env%sab_core)
CALL release_neighbor_list_sets(ks_env%sab_xb)
CALL release_neighbor_list_sets(ks_env%sab_xtb_nonbond)
CALL release_neighbor_list_sets(ks_env%sab_all)
CALL release_neighbor_list_sets(ks_env%sab_aux_fit)
CALL release_neighbor_list_sets(ks_env%sab_aux_fit_vs_orb)

52
src/qs_neighbor_lists.F
View file

@ -87,7 +87,8 @@ MODULE qs_neighbor_lists
get_iterator_task, neighbor_list_iterate, neighbor_list_iterator_create, &
neighbor_list_iterator_p_type, neighbor_list_iterator_release, neighbor_list_p_type, &
neighbor_list_set_p_type, neighbor_list_set_type, release_neighbor_list_sets
USE string_utilities, ONLY: compress
USE string_utilities, ONLY: compress,&
uppercase
USE subcell_types, ONLY: allocate_subcell,&
deallocate_subcell,&
give_ijk_subcell,&
@ -295,15 +296,16 @@ CONTAINS
CHARACTER(len=*), PARAMETER :: routineN = 'build_qs_neighbor_lists'
CHARACTER(LEN=2) :: element_symbol, element_symbol2
CHARACTER(LEN=default_string_length) :: print_key_path
INTEGER :: handle, ikind, iw, maxatom, nkind, &
numshells, zat
INTEGER :: handle, ikind, ingp, iw, jkind, maxatom, &
ngp, nkind, numshells, zat
LOGICAL :: all_potential_present, almo, dftb, do_hfx, dokp, explicit, gth_potential_present, &
lri_optbas, lrigpw, mic, molecule_only, nddo, paw_atom, paw_atom_present, rigpw, &
sgp_potential_present, xtb
LOGICAL, ALLOCATABLE, DIMENSION(:) :: all_present, aux_fit_present, aux_present, &
core_present, default_present, oce_present, orb_present, ppl_present, ppnl_present, &
ri_present, xb1_atom, xb2_atom
core_present, default_present, nonbond1_atom, nonbond2_atom, oce_present, orb_present, &
ppl_present, ppnl_present, ri_present, xb1_atom, xb2_atom
REAL(dp) :: almo_rcov, almo_rvdw, rcut, roperator, &
subcells
REAL(dp), ALLOCATABLE, DIMENSION(:) :: all_pot_rad, aux_fit_radius, c_radius, calpha, &
@ -325,8 +327,8 @@ CONTAINS
TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set
TYPE(neighbor_list_set_p_type), DIMENSION(:), POINTER :: saa_list, sab_all, sab_almo, &
sab_aux_fit, sab_aux_fit_asymm, sab_aux_fit_vs_orb, sab_cn, sab_core, sab_gcp, sab_kp, &
sab_lrc, sab_orb, sab_scp, sab_se, sab_tbe, sab_vdw, sab_xb, sac_ae, sac_lri, sac_ppl, &
sap_oce, sap_ppnl, soa_list, soo_list
sab_lrc, sab_orb, sab_scp, sab_se, sab_tbe, sab_vdw, sab_xb, sab_xtb_nonbond, sac_ae, &
sac_lri, sac_ppl, sap_oce, sap_ppnl, soa_list, soo_list
TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
TYPE(paw_proj_set_type), POINTER :: paw_proj
TYPE(qs_dftb_atom_type), POINTER :: dftb_atom
@ -356,6 +358,7 @@ CONTAINS
NULLIFY (sab_tbe)
NULLIFY (sab_core)
NULLIFY (sab_xb)
NULLIFY (sab_xtb_nonbond)
NULLIFY (sab_all)
NULLIFY (sab_vdw)
NULLIFY (sab_cn)
@ -401,6 +404,7 @@ CONTAINS
sab_tbe=sab_tbe, &
sab_core=sab_core, &
sab_xb=sab_xb, &
sab_xtb_nonbond=sab_xtb_nonbond, &
sab_scp=sab_scp, &
sab_all=sab_all, &
sab_almo=sab_almo, &
@ -801,6 +805,40 @@ CONTAINS
CALL build_neighbor_lists(sab_xb, particle_set, atom2d, cell, pair_radius, &
symmetric=.FALSE., subcells=subcells, operator_type="PP", nlname="sab_xb")
CALL set_ks_env(ks_env=ks_env, sab_xb=sab_xb)
CALL write_neighbor_lists(sab_xb, particle_set, cell, para_env, neighbor_list_section, &
"/SAB_XB", "sab_xb", "XB bonding")
! nonbonded interactions list
IF (dft_control%qs_control%xtb_control%do_nonbonded) THEN
ngp = SIZE(dft_control%qs_control%xtb_control%nonbonded%pot)
ALLOCATE (nonbond1_atom(nkind), nonbond2_atom(nkind))
nonbond1_atom = .FALSE.
nonbond2_atom = .FALSE.
DO ingp = 1, ngp
DO ikind = 1, nkind
rcut = SQRT(dft_control%qs_control%xtb_control%nonbonded%pot(ingp)%pot%rcutsq)
c_radius = rcut
CALL get_atomic_kind(atomic_kind_set(ikind), element_symbol=element_symbol)
CALL uppercase(element_symbol)
IF (TRIM(dft_control%qs_control%xtb_control%nonbonded%pot(ingp)%pot%at1) == TRIM(element_symbol)) THEN
nonbond1_atom(ikind) = .TRUE.
DO jkind = 1, nkind
CALL get_atomic_kind(atomic_kind_set(jkind), element_symbol=element_symbol2)
CALL uppercase(element_symbol2)
IF (TRIM(dft_control%qs_control%xtb_control%nonbonded%pot(ingp)%pot%at2) == TRIM(element_symbol2)) THEN
nonbond2_atom(jkind) = .TRUE.
ENDIF
ENDDO
END IF
END DO
CALL pair_radius_setup(nonbond1_atom, nonbond2_atom, c_radius, c_radius, pair_radius)
CALL build_neighbor_lists(sab_xtb_nonbond, particle_set, atom2d, cell, pair_radius, &
symmetric=.FALSE., subcells=subcells, operator_type="PP", nlname="sab_xtb_nonbond")
CALL set_ks_env(ks_env=ks_env, sab_xtb_nonbond=sab_xtb_nonbond)
CALL write_neighbor_lists(sab_xtb_nonbond, particle_set, cell, para_env, neighbor_list_section, &
"/SAB_XTB_NONBOND", "sab_xtb_nonbond", "XTB NONBONDED INTERACTIONS")
END DO
ENDIF
END IF
! Build the neighbor lists for the vdW pair potential

6
src/qs_scf_output.F
View file

@ -372,6 +372,12 @@ CONTAINS
"Electronic energy: ", energy%hartree, &
"DFTB3 3rd order energy: ", energy%dftb3, &
"Dispersion energy: ", energy%dispersion
IF (dft_control%qs_control%xtb_control%xb_interaction) &
WRITE (UNIT=output_unit, FMT="(T3,A,T56,F25.14)") &
"Correction for halogen bonding: ", energy%xtb_xb_inter
IF (dft_control%qs_control%xtb_control%do_nonbonded) &
WRITE (UNIT=output_unit, FMT="(T3,A,T56,F25.14)") &
"Correction for nonbonded interactions: ", energy%xtb_nonbonded
IF (energy%efield /= 0.0_dp) &
WRITE (UNIT=output_unit, FMT="(T3,A,T56,F25.14)") &
"Electric field interaction energy: ", energy%efield

289
src/xtb_matrices.F
View file

@ -32,7 +32,8 @@ MODULE xtb_matrices
USE cp_dbcsr_output, ONLY: cp_dbcsr_write_sparse_matrix
USE cp_log_handling, ONLY: cp_get_default_logger,&
cp_logger_get_default_io_unit,&
cp_logger_type
cp_logger_type,&
cp_to_string
USE cp_output_handling, ONLY: cp_p_file,&
cp_print_key_finished_output,&
cp_print_key_should_output,&
@ -42,16 +43,31 @@ MODULE xtb_matrices
dbcsr_add, dbcsr_copy, dbcsr_create, dbcsr_dot, dbcsr_finalize, dbcsr_get_block_p, &
dbcsr_multiply, dbcsr_p_type, dbcsr_type
USE efield_tb_methods, ONLY: efield_tb_matrix
USE ewald_environment_types, ONLY: ewald_env_get,&
ewald_environment_type
USE fparser, ONLY: evalfd,&
finalizef
USE input_section_types, ONLY: section_vals_get_subs_vals,&
section_vals_type,&
section_vals_val_get
USE kinds, ONLY: dp
USE kinds, ONLY: default_string_length,&
dp
USE kpoint_types, ONLY: get_kpoint_info,&
kpoint_type
USE message_passing, ONLY: mp_minloc,&
mp_sum
USE mulliken, ONLY: ao_charges
USE orbital_pointers, ONLY: ncoset
USE pair_potential, ONLY: init_genpot
USE pair_potential_types, ONLY: not_initialized,&
pair_potential_p_type,&
pair_potential_pp_create,&
pair_potential_pp_release,&
pair_potential_pp_type,&
pair_potential_single_clean,&
pair_potential_single_copy,&
pair_potential_single_type
USE pair_potential_util, ONLY: ener_pot
USE particle_types, ONLY: particle_type
USE qs_charge_mixing, ONLY: charge_mixing
USE qs_condnum, ONLY: overlap_condnum
@ -83,6 +99,8 @@ MODULE xtb_matrices
USE qs_rho_types, ONLY: qs_rho_get,&
qs_rho_type
USE qs_scf_types, ONLY: qs_scf_env_type
USE string_utilities, ONLY: compress,&
uppercase
USE virial_methods, ONLY: virial_pair_force
USE virial_types, ONLY: virial_type
USE xtb_coulomb, ONLY: build_xtb_coulomb
@ -132,14 +150,13 @@ CONTAINS
npgfb, nsgfa, nsgfb
INTEGER, DIMENSION(:, :), POINTER :: first_sgfa, first_sgfb
INTEGER, DIMENSION(:, :, :), POINTER :: cell_to_index
LOGICAL :: defined, diagblock, floating_a, found, &
ghost_a, norml1, norml2, omit_headers, &
use_arnoldi, use_virial, xb_inter
LOGICAL :: defined, diagblock, do_nonbonded, floating_a, found, ghost_a, norml1, norml2, &
omit_headers, use_arnoldi, use_virial, xb_inter
LOGICAL, ALLOCATABLE, DIMENSION(:) :: floating, ghost
REAL(KIND=dp) :: alphaa, alphab, derepij, dfp, dhij, dr, drk, drx, ena, enb, erep, erepij, &
etaa, etab, exb, f0, f1, fen, fhua, fhub, fhud, foab, hij, k2sh, kab, kcnd, kcnp, kcns, &
kd, ken, kf, kg, kia, kjb, kp, ks, ksp, kx2, kxr, rcova, rcovab, rcovb, rrab, zneffa, &
zneffb
REAL(KIND=dp) :: alphaa, alphab, derepij, dfp, dhij, dr, drk, drx, ena, enb, enonbonded, &
erep, erepij, etaa, etab, exb, f0, f1, fen, fhua, fhub, fhud, foab, hij, k2sh, kab, kcnd, &
kcnp, kcns, kd, ken, kf, kg, kia, kjb, kp, ks, ksp, kx2, kxr, rcova, rcovab, rcovb, rrab, &
zneffa, zneffb
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: cnumbers, kx
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: dfblock, huckel, kcnlk, owork, rcab
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :) :: oint, sint
@ -168,7 +185,7 @@ CONTAINS
TYPE(neighbor_list_iterator_p_type), &
DIMENSION(:), POINTER :: nl_iterator
TYPE(neighbor_list_set_p_type), DIMENSION(:), &
POINTER :: sab_orb, sab_xb
POINTER :: sab_orb, sab_xb, sab_xtb_nonbond
TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
TYPE(qs_dispersion_type), POINTER :: dispersion_env
TYPE(qs_energy_type), POINTER :: energy
@ -180,6 +197,8 @@ CONTAINS
TYPE(xtb_atom_type), POINTER :: xtb_atom_a, xtb_atom_b
TYPE(xtb_control_type), POINTER :: xtb_control
!, na1
CALL timeset(routineN, handle)
NULLIFY (logger, virial, atprop)
@ -202,12 +221,16 @@ CONTAINS
nkind = SIZE(atomic_kind_set)
xtb_control => dft_control%qs_control%xtb_control
xb_inter = xtb_control%xb_interaction
do_nonbonded = xtb_control%do_nonbonded
nimg = dft_control%nimages
nderivatives = 0
IF (calculate_forces) nderivatives = 1
IF (dft_control%tddfpt2_control%enabled) nderivatives = 1
maxder = ncoset(nderivatives)
IF (xb_inter) energy%xtb_xb_inter = 0.0_dp
IF (do_nonbonded) energy%xtb_nonbonded = 0.0_dp
! global parameters (Table 2 from Ref.)
ks = xtb_control%ks
kp = xtb_control%kp
@ -842,8 +865,19 @@ CONTAINS
CALL xb_interaction(exb, neighbor_atoms, atom_of_kind, kind_of, sab_xb, kx, kx2, rcab, &
calculate_forces, use_virial, force, virial, atprop)
END IF
enonbonded = 0.0_dp
IF (do_nonbonded) THEN
! nonbonded interactions
NULLIFY (sab_xtb_nonbond)
CALL get_qs_env(qs_env=qs_env, sab_xtb_nonbond=sab_xtb_nonbond)
CALL nonbonded_correction(enonbonded, force, qs_env, xtb_control, sab_xtb_nonbond, &
atomic_kind_set, calculate_forces, use_virial, virial, atprop, atom_of_kind)
ENDIF
! set repulsive energy
erep = erep + exb
erep = erep + exb + enonbonded
IF (xb_inter) energy%xtb_xb_inter = exb
IF (do_nonbonded) energy%xtb_nonbonded = enonbonded
CALL mp_sum(erep, para_env%group)
energy%repulsive = erep
@ -1148,6 +1182,12 @@ CONTAINS
"Zeroth order Hamiltonian energy: ", energy%core, &
"Charge fluctuation energy: ", energy%hartree, &
"London dispersion energy: ", energy%dispersion
IF (dft_control%qs_control%xtb_control%xb_interaction) &
WRITE (UNIT=iounit, FMT="(T9,A,T60,F20.10)") &
"Correction for halogen bonding: ", energy%xtb_xb_inter
IF (dft_control%qs_control%xtb_control%do_nonbonded) &
WRITE (UNIT=iounit, FMT="(T9,A,T60,F20.10)") &
"Correction for nonbonded interactions: ", energy%xtb_nonbonded
IF (ABS(energy%efield) > 1.e-12_dp) THEN
WRITE (UNIT=iounit, FMT="(T9,A,T60,F20.10)") &
"Electric field interaction energy: ", energy%efield
@ -1228,6 +1268,233 @@ CONTAINS
END SUBROUTINE xb_neighbors
! **************************************************************************************************
!> \brief Computes a correction for nonbonded interactions based on a generic potential
!>
!> \param enonbonded energy contribution
!> \param force ...
!> \param qs_env ...
!> \param xtb_control ...
!> \param sab_xtb_nonbond ...
!> \param atomic_kind_set ...
!> \param calculate_forces ...
!> \param use_virial ...
!> \param virial ...
!> \param atprop ...
!> \param atom_of_kind ..
!> \par History
!> 12.2018 JGH
!> \version 1.1
! **************************************************************************************************
SUBROUTINE nonbonded_correction(enonbonded, force, qs_env, xtb_control, sab_xtb_nonbond, &
atomic_kind_set, calculate_forces, use_virial, virial, atprop, atom_of_kind)
REAL(dp), INTENT(INOUT) :: enonbonded
TYPE(qs_force_type), DIMENSION(:), INTENT(INOUT), &
POINTER :: force
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(xtb_control_type), POINTER :: xtb_control
TYPE(neighbor_list_set_p_type), DIMENSION(:), &
INTENT(IN), POINTER :: sab_xtb_nonbond
TYPE(atomic_kind_type), DIMENSION(:), INTENT(IN), &
POINTER :: atomic_kind_set
LOGICAL, INTENT(IN) :: calculate_forces, use_virial
TYPE(virial_type), INTENT(IN), POINTER :: virial
TYPE(atprop_type), INTENT(IN), POINTER :: atprop
INTEGER, DIMENSION(:), INTENT(IN) :: atom_of_kind
CHARACTER(len=*), PARAMETER :: routineN = 'nonbonded_correction'
CHARACTER(LEN=default_string_length) :: def_error, this_error
INTEGER :: atom_i, atom_j, handle, iatom, ikind, &
jatom, jkind, kk, ntype
INTEGER, DIMENSION(3) :: cell
LOGICAL :: do_ewald
REAL(KIND=dp) :: dedf, dr, dx, energy_cutoff, err, lerr, &
rcut
REAL(KIND=dp), DIMENSION(3) :: fij, rij
TYPE(ewald_environment_type), POINTER :: ewald_env
TYPE(neighbor_list_iterator_p_type), &
DIMENSION(:), POINTER :: nl_iterator
TYPE(pair_potential_p_type), POINTER :: nonbonded
TYPE(pair_potential_pp_type), POINTER :: potparm
TYPE(pair_potential_single_type), POINTER :: pot
TYPE(section_vals_type), POINTER :: nonbonded_section
!
CALL timeset(routineN, handle)
NULLIFY (nonbonded)
NULLIFY (potparm)
NULLIFY (ewald_env)
nonbonded => xtb_control%nonbonded
do_ewald = xtb_control%do_ewald
CALL get_qs_env(qs_env=qs_env, ewald_env=ewald_env)
ntype = SIZE(atomic_kind_set)
CALL pair_potential_pp_create(potparm, ntype)
!Assign input and potential info to potparm_nonbond
CALL force_field_pack_nonbond_pot_correction(atomic_kind_set, nonbonded, potparm, ewald_env, do_ewald)
!Initialize genetic potential
CALL init_genpot(potparm, ntype)
NULLIFY (pot)
enonbonded = 0._dp
energy_cutoff = 0._dp
CALL neighbor_list_iterator_create(nl_iterator, sab_xtb_nonbond)
DO WHILE (neighbor_list_iterate(nl_iterator) == 0)
CALL get_iterator_info(nl_iterator, ikind=ikind, jkind=jkind, &
iatom=iatom, jatom=jatom, r=rij, cell=cell)
pot => potparm%pot(ikind, jkind)%pot
dr = SQRT(rij(1)**2 + rij(2)**2 + rij(3)**2)
rcut = SQRT(pot%rcutsq)
IF (dr <= rcut) THEN
! splines not implemented
enonbonded = enonbonded + ener_pot(pot, dr, energy_cutoff)
IF (atprop%energy) THEN
atprop%atecc(iatom) = atprop%atecc(iatom) + 0.5_dp*ener_pot(pot, dr, energy_cutoff)
atprop%atecc(jatom) = atprop%atecc(jatom) + 0.5_dp*ener_pot(pot, dr, energy_cutoff)
END IF
ENDIF
IF (calculate_forces) THEN
kk = SIZE(pot%type)
IF (kk /= 1) THEN
CALL cp_warn(__LOCATION__, "Generic potential with type > 1 not implemented.")
CPABORT("pot type")
ENDIF
! rmin and rmax and rcut
IF ((pot%set(kk)%rmin /= not_initialized) .AND. (dr < pot%set(kk)%rmin)) CYCLE
! An upper boundary for the potential definition was defined
IF ((pot%set(kk)%rmax /= not_initialized) .AND. (dr >= pot%set(kk)%rmax)) CYCLE
! If within limits let's compute the potential...
IF (dr <= rcut) THEN
NULLIFY (nonbonded_section)
nonbonded_section => section_vals_get_subs_vals(qs_env%input, "DFT%QS%xTB%NONBONDED")
CALL section_vals_val_get(nonbonded_section, "DX", r_val=dx)
CALL section_vals_val_get(nonbonded_section, "ERROR_LIMIT", r_val=lerr)
dedf = evalfd(pot%set(kk)%gp%myid, 1, pot%set(kk)%gp%values, dx, err)
IF (ABS(err) > lerr) THEN
WRITE (this_error, "(A,G12.6,A)") "(", err, ")"
WRITE (def_error, "(A,G12.6,A)") "(", lerr, ")"
CALL compress(this_error, .TRUE.)
CALL compress(def_error, .TRUE.)
CALL cp_warn(__LOCATION__, &
'ASSERTION (cond) failed at line '//cp_to_string(__LINE__)// &
' Error '//TRIM(this_error)//' in computing numerical derivatives larger then'// &
TRIM(def_error)//' .')
END IF
!
atom_i = atom_of_kind(iatom)
atom_j = atom_of_kind(jatom)
!
fij(1:3) = dedf*rij(1:3)/pot%set(kk)%gp%values(1)
force(ikind)%repulsive(:, atom_i) = force(ikind)%repulsive(:, atom_i) - fij(:)
force(jkind)%repulsive(:, atom_j) = force(jkind)%repulsive(:, atom_j) + fij(:)
IF (use_virial) THEN
CALL virial_pair_force(virial%pv_virial, -1._dp, fij, rij)
IF (atprop%stress) THEN
CALL virial_pair_force(atprop%atstress(:, :, iatom), -0.5_dp, fij, rij)
CALL virial_pair_force(atprop%atstress(:, :, jatom), -0.5_dp, fij, rij)
END IF
END IF
END IF
END IF
NULLIFY (pot)
END DO
CALL neighbor_list_iterator_release(nl_iterator)
CALL finalizef()
IF (ASSOCIATED(potparm)) THEN
CALL pair_potential_pp_release(potparm)
ENDIF
CALL timestop(handle)
!
END SUBROUTINE nonbonded_correction
! **************************************************************************************************
!> \brief ...
!> \param atomic_kind_set ...
!> \param nonbonded ...
!> \param potparm ...
!> \param ewald_env ...
!> \param do_ewald ...
! **************************************************************************************************
SUBROUTINE force_field_pack_nonbond_pot_correction(atomic_kind_set, nonbonded, potparm, ewald_env, do_ewald)
! routine based on force_field_pack_nonbond
TYPE(atomic_kind_type), DIMENSION(:), INTENT(IN), &
POINTER :: atomic_kind_set
TYPE(pair_potential_p_type), INTENT(IN), POINTER :: nonbonded
TYPE(pair_potential_pp_type), INTENT(INOUT), &
POINTER :: potparm
TYPE(ewald_environment_type), INTENT(IN), POINTER :: ewald_env
LOGICAL, INTENT(IN) :: do_ewald
CHARACTER(LEN=default_string_length) :: name_atm_a, name_atm_a_local, &
name_atm_b, name_atm_b_local
INTEGER :: ikind, ingp, iw, jkind
LOGICAL :: found
REAL(KIND=dp) :: ewald_rcut
TYPE(atomic_kind_type), POINTER :: atomic_kind
TYPE(pair_potential_single_type), POINTER :: pot
NULLIFY (pot)
DO ikind = 1, SIZE(atomic_kind_set)
atomic_kind => atomic_kind_set(ikind)
CALL get_atomic_kind(atomic_kind=atomic_kind, name=name_atm_a_local)
DO jkind = ikind, SIZE(atomic_kind_set)
atomic_kind => atomic_kind_set(jkind)
CALL get_atomic_kind(atomic_kind=atomic_kind, name=name_atm_b_local)
found = .FALSE.
name_atm_a = name_atm_a_local
name_atm_b = name_atm_b_local
CALL uppercase(name_atm_a)
CALL uppercase(name_atm_b)
pot => potparm%pot(ikind, jkind)%pot
IF (ASSOCIATED(nonbonded)) THEN
DO ingp = 1, SIZE(nonbonded%pot)
IF ((TRIM(nonbonded%pot(ingp)%pot%at1) == "*") .OR. &
(TRIM(nonbonded%pot(ingp)%pot%at2) == "*")) CYCLE
IF (iw > 0) WRITE (iw, *) "TESTING ", TRIM(name_atm_a), TRIM(name_atm_b), &
" with ", TRIM(nonbonded%pot(ingp)%pot%at1), &
TRIM(nonbonded%pot(ingp)%pot%at2)
IF ((((name_atm_a) == (nonbonded%pot(ingp)%pot%at1)) .AND. &
((name_atm_b) == (nonbonded%pot(ingp)%pot%at2))) .OR. &
(((name_atm_b) == (nonbonded%pot(ingp)%pot%at1)) .AND. &
((name_atm_a) == (nonbonded%pot(ingp)%pot%at2)))) THEN
CALL pair_potential_single_copy(nonbonded%pot(ingp)%pot, pot)
! multiple potential not implemented, simply overwriting
IF (found) &
CALL cp_warn(__LOCATION__, &
"Multiple NONBONDED declaration: "//TRIM(name_atm_a)// &
" and "//TRIM(name_atm_b)//" OVERWRITING! ")
END IF
IF (iw > 0) WRITE (iw, *) " FOUND ", TRIM(name_atm_a), " ", TRIM(name_atm_b)
found = .TRUE.
END DO
ENDIF
IF (.NOT. found) THEN
CALL pair_potential_single_clean(pot)
ENDIF
ENDDO !jkind
ENDDO !ikind
! Cutoff is defined always as the maximum between the FF and Ewald
IF (do_ewald) THEN
CALL ewald_env_get(ewald_env, rcut=ewald_rcut)
pot%rcutsq = MAX(pot%rcutsq, ewald_rcut*ewald_rcut)
ENDIF
END SUBROUTINE force_field_pack_nonbond_pot_correction
! **************************************************************************************************
!> \brief Computes the interaction term between Br/I/At and donor atoms
!>
!> \param exb ...

3
tests/xTB/regtest-debug/TEST_FILES
View file

@ -9,4 +9,7 @@ ch2o_t03.inp 86 1.0E-12 0.11459198
ch2o_t04.inp 86 1.0E-12 0.109829662200E+01
ch2o_t05.inp 87 1.0E-12 0.226843652140E+02
ch2o_t06.inp 87 1.0E-12 0.230117849452E+02
ch3br_nonbond_2.inp 1 1.0E-12 -11.78062819229165
ch3br_nonbond.inp 1 1.0E-12 -11.78062802267188
ch3br_atprop_nonbond.inp 1 1.0E-10 -11.78101985886981
#EOF

71
tests/xTB/regtest-debug/ch3br_atprop_nonbond.inp Normal file
View file

@ -0,0 +1,71 @@
&FORCE_EVAL
STRESS_TENSOR ANALYTICAL
&PROPERTIES
&ATOMIC
ENERGY
PRESSURE
&END ATOMIC
&END PROPERTIES
METHOD Quickstep
&DFT
&QS
METHOD xTB
EPS_DEFAULT 1.0E-12
&xTB
DO_EWALD T
CHECK_ATOMIC_CHARGES OFF
DO_NONBONDED T
&NONBONDED
&GENPOT
ATOMS Br Xe
FUNCTION Aparam*exp(-Bparam*r)-Cparam/r**8
PARAMETERS Aparam Bparam Cparam
VALUES 49.507 1.441 6703.192
VARIABLES r
&END GENPOT
DX 0.0001
ERROR_LIMIT 0.0000001
&END NONBONDED
&END XTB
&END QS
&SCF
SCF_GUESS RESTART
EPS_SCF 1.0E-8
MAX_SCF 50
&OT
MINIMIZER DIIS
PRECONDITIONER FULL_S_INVERSE
&END
&OUTER_SCF
EPS_SCF 1.0E-8
MAX_SCF 10
&END
&END SCF
&END DFT
&SUBSYS
&CELL
ABC 20.0 20.0 20.0
&END CELL
&COORD
C -0.000001 -1.540826 -0.000000
H 1.041396 -1.871752 -0.000000
H -0.520681 -1.871772 0.901873
H -0.520681 -1.871772 -0.901873
Br -0.000001 0.424579 0.000000
Xe -0.000000 3.0 -0.000000
&END COORD
&END SUBSYS
&END FORCE_EVAL
&GLOBAL
PROJECT CH3Br-Xe
PRINT_LEVEL LOW
RUN_TYPE MD
&END GLOBAL
&MOTION
&MD
ENSEMBLE NVE
STEPS 3
TIMESTEP 0.48
TEMPERATURE 330.0
&END MD
&END MOTION

56
tests/xTB/regtest-debug/ch3br_nonbond.inp Normal file
View file

@ -0,0 +1,56 @@
&FORCE_EVAL
METHOD Quickstep
&DFT
&QS
METHOD xTB
EPS_DEFAULT 1.0E-10
EPS_PGF_ORB 1.0E-20
&xTB
DO_EWALD F
CHECK_ATOMIC_CHARGES OFF
DO_NONBONDED T
&NONBONDED
&GENPOT
ATOMS Br Xe
FUNCTION Aparam*exp(-Bparam*r)-Cparam/r**8
PARAMETERS Aparam Bparam Cparam
VALUES 49.507 1.441 6703.192
VARIABLES r
&END GENPOT
&END NONBONDED
&END XTB
&END QS
&SCF
SCF_GUESS RESTART
EPS_SCF 1.0E-5
MAX_SCF 100
&OT
MINIMIZER DIIS
PRECONDITIONER FULL_S_INVERSE
&END
&OUTER_SCF
EPS_SCF 1.0E-5
MAX_SCF 900
&END
&END SCF
&END DFT
&SUBSYS
&CELL
ABC 10 10 10
PERIODIC NONE
&END CELL
&COORD
C -0.000001 -1.540826 -0.000000
H 1.041396 -1.871752 -0.000000
H -0.520681 -1.871772 0.901873
H -0.520681 -1.871772 -0.901873
Br -0.000001 0.424579 0.000000
Xe -0.000000 3.0 -0.000000
&END COORD
&END SUBSYS
&END FORCE_EVAL
&GLOBAL
PROJECT CH3Br-Xe
PRINT_LEVEL LOW
RUN_TYPE ENERGY
&END GLOBAL

70
tests/xTB/regtest-debug/ch3br_nonbond_2.inp Normal file
View file

@ -0,0 +1,70 @@
&FORCE_EVAL
METHOD Quickstep
&DFT
&QS
METHOD xTB
&xTB
DO_EWALD F
CHECK_ATOMIC_CHARGES OFF
DO_NONBONDED T
&NONBONDED
&GENPOT
ATOMS Br Xe
FUNCTION Aparam*exp(-Bparam*r)-Cparam/r**8
PARAMETERS Aparam Bparam Cparam
VALUES 49.507 1.441 6703.192
VARIABLES r
&END GENPOT
&END NONBONDED
&END XTB
&END QS
&SCF
SCF_GUESS RESTART
EPS_SCF 1.0E-5
MAX_SCF 100
&OT
MINIMIZER DIIS
PRECONDITIONER FULL_S_INVERSE
&END
&OUTER_SCF
EPS_SCF 1.0E-5
MAX_SCF 900
&END
&END SCF
&EFIELD
&END
&PRINT
&MOMENTS ON
PERIODIC .FALSE.
REFERENCE COM
&END
&END
&END DFT
&SUBSYS
&CELL
ABC 10 10 10
PERIODIC NONE
&END CELL
&COORD
C -0.000001 -1.540826 -0.000000
H 1.041396 -1.871752 -0.000000
H -0.520681 -1.871772 0.901873
H -0.520681 -1.871772 -0.901873
Br -0.000001 0.424579 0.000000
Xe -0.000000 3.0 -0.000000
&END COORD
&END SUBSYS
&END FORCE_EVAL
&GLOBAL
PROJECT CH3Br-Xe
PRINT_LEVEL LOW
RUN_TYPE DEBUG
&END GLOBAL
&DEBUG
DEBUG_FORCES F
DEBUG_STRESS_TENSOR F
DEBUG_DIPOLE T
DEBUG_POLARIZABILITY F
STOP_ON_MISMATCH T
DE 0.0002
&END

176
tests/xTB/regtest-debug/h2o-atprop_nonbonded.inp Normal file
View file

@ -0,0 +1,176 @@
&FORCE_EVAL
STRESS_TENSOR ANALYTICAL
&PROPERTIES
&ATOMIC
ENERGY
PRESSURE
&END ATOMIC
&END PROPERTIES
&DFT
&QS
METHOD xTB
&xTB
DO_EWALD T
COULOMB_INTERACTION T
TB3_INTERACTION T
DO_NONBONDED T
&NONBONDED
&GENPOT
ATOMS O O
FUNCTION 4*epsilon*((sigma/(r+sigma))**12 - (sigma/(r+sigma))**6)
PARAMETERS epsilon sigma
VALUES 0.00380879917604 2.31204899506
VARIABLES r
RCUT 4.0
RMIN 1.0
RMAX 4.0
&END GENPOT
DX 0.001
ERROR_LIMIT 0.0000000001
&END NONBONDED
&END
&END QS
&SCF
SCF_GUESS MOPAC
EPS_SCF 1.E-5
&MIXING
METHOD DIRECT_P_MIXING
ALPHA 0.75
&END
MAX_SCF 20
&END SCF
&POISSON
&EWALD
O_SPLINE 6
&END EWALD
&END POISSON
&END DFT
&SUBSYS
&CELL
ABC 9.865 9.865 9.865
&END CELL
&COORD
O -4.583 5.333 1.560 H2O
H -3.777 5.331 0.943 H2O
H -5.081 4.589 1.176 H2O
O -0.083 4.218 0.070 H2O
H -0.431 3.397 0.609 H2O
H 0.377 3.756 -0.688 H2O
O -1.488 2.692 4.125 H2O
H -2.465 2.433 3.916 H2O
H -1.268 2.145 4.952 H2O
O -2.461 -2.548 -6.136 H2O
H -1.892 -2.241 -6.921 H2O
H -1.970 -3.321 -5.773 H2O
O 4.032 0.161 2.183 H2O
H 4.272 -0.052 1.232 H2O
H 4.044 -0.760 2.641 H2O
O 2.950 -3.497 -1.006 H2O
H 2.599 -3.901 -0.129 H2O
H 3.193 -4.283 -1.533 H2O
O -2.890 -4.797 -2.735 H2O
H -2.810 -5.706 -2.297 H2O
H -2.437 -4.128 -2.039 H2O
O -0.553 0.922 -3.731 H2O
H -0.163 1.552 -3.085 H2O
H -1.376 0.544 -3.469 H2O
O 4.179 4.017 4.278 H2O
H 3.275 3.832 3.876 H2O
H 4.658 4.492 3.572 H2O
O 5.739 1.425 3.944 H2O
H 5.125 2.066 4.346 H2O
H 5.173 1.181 3.097 H2O
O 0.988 -0.760 -5.445 H2O
H 1.640 -1.372 -4.989 H2O
H 0.546 -0.220 -4.762 H2O
O -0.748 1.985 1.249 H2O
H -0.001 1.490 1.540 H2O
H -1.160 2.255 2.109 H2O
O 4.127 -0.234 -3.149 H2O
H 5.022 -0.436 -3.428 H2O
H 3.540 -0.918 -3.601 H2O
O -2.473 2.768 -1.395 H2O
H -1.533 2.719 -1.214 H2O
H -2.702 1.808 -1.479 H2O
O -0.124 -2.116 2.404 H2O
H 0.612 -2.593 2.010 H2O
H 0.265 -1.498 3.089 H2O
O 0.728 2.823 -2.190 H2O
H 0.646 3.694 -2.685 H2O
H 1.688 2.705 -1.947 H2O
O 4.256 -5.427 -2.644 H2O
H 5.222 -5.046 -2.479 H2O
H 4.174 -5.628 -3.593 H2O
O -3.178 -0.508 -4.227 H2O
H -2.762 -1.221 -4.818 H2O
H -3.603 0.073 -4.956 H2O
O -1.449 5.300 -4.805 H2O
H -1.397 4.470 -5.317 H2O
H -2.102 5.091 -4.067 H2O
O 3.354 2.192 -1.755 H2O
H 3.407 1.433 -2.405 H2O
H 3.971 2.958 -2.196 H2O
O 1.773 -4.018 1.769 H2O
H 1.121 -4.532 1.201 H2O
H 1.975 -4.529 2.618 H2O
O 1.526 1.384 2.712 H2O
H 2.317 1.070 2.251 H2O
H 1.353 0.657 3.364 H2O
O 2.711 -2.398 -4.253 H2O
H 2.202 -3.257 -4.120 H2O
H 3.305 -2.610 -5.099 H2O
O 6.933 0.093 -1.393 H2O
H 6.160 -0.137 -0.795 H2O
H 6.748 -0.394 -2.229 H2O
O -5.605 -2.549 3.151 H2O
H -4.756 -2.503 3.616 H2O
H -5.473 -3.187 2.378 H2O
O 0.821 -4.406 6.516 H2O
H 0.847 -3.675 7.225 H2O
H -0.014 -4.240 5.988 H2O
O 1.577 3.933 3.762 H2O
H 1.221 2.975 3.640 H2O
H 1.367 4.126 4.659 H2O
O -2.111 -3.741 -0.219 H2O
H -1.378 -4.425 -0.036 H2O
H -1.825 -2.775 0.003 H2O
O 0.926 -1.961 -2.063 H2O
H 0.149 -1.821 -1.402 H2O
H 1.725 -2.303 -1.536 H2O
O 4.531 -1.030 -0.547 H2O
H 4.290 -1.980 -0.581 H2O
H 4.292 -0.597 -1.390 H2O
O -0.740 -1.262 -0.029 H2O
H -1.272 -0.422 -0.099 H2O
H -0.403 -1.349 0.873 H2O
O 3.655 3.021 0.988 H2O
H 2.706 3.053 1.282 H2O
H 3.542 2.615 0.020 H2O
&END COORD
&END SUBSYS
&END FORCE_EVAL
&GLOBAL
PROJECT h2o
RUN_TYPE MD
PRINT_LEVEL LOW
&END GLOBAL
&MOTION
&MD
ENSEMBLE NPT_I
STEPS 3
TIMESTEP 0.48
TEMPERATURE 330.0
&THERMOSTAT
&NOSE
LENGTH 3
YOSHIDA 3
TIMECON [wavenumber_t] 2000
MTS 2
&END NOSE
&END
&BAROSTAT
PRESSURE 1.0
TIMECON 300
&END BAROSTAT
&END MD
&END MOTION