Skip to content

Commit

Permalink
Write xTB parameters to output
Browse files Browse the repository at this point in the history
  • Loading branch information
@juerghutter
juerghutter committed Mar 19, 2019
1 parent bac399e commit e5cf0c8
Show file tree
Hide file tree
Showing 4 changed files with 153 additions and 15 deletions.
80 changes: 72 additions & 8 deletions src/cp_control_utils.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -7,16 +7,15 @@
!> \brief Utilities to set up the control types
! **************************************************************************************************
MODULE cp_control_utils
!qxtb citation
USE bibliography, ONLY: &
Andreussi2012, Dewar1977, Dewar1985, Elstner1998, Fattebert2002, Hu2007, Krack2000, &
Lippert1997, Lippert1999, Porezag1995, Repasky2002, Rocha2006, Schenter2008, Seifert1996, &
Souza2002, Stengel2009, Stewart1989, Stewart2007, Thiel1992, Umari2002, VanVoorhis2015, &
VandeVondele2005a, VandeVondele2005b, Zhechkov2005, cite_reference
Andreussi2012, Dewar1977, Dewar1985, Elstner1998, Fattebert2002, Grimme2017, Hu2007, &
Krack2000, Lippert1997, Lippert1999, Porezag1995, Repasky2002, Rocha2006, Schenter2008, &
Seifert1996, Souza2002, Stengel2009, Stewart1989, Stewart2007, Thiel1992, Umari2002, &
VanVoorhis2015, VandeVondele2005a, VandeVondele2005b, Zhechkov2005, cite_reference
USE cp_control_types, ONLY: &
admm_control_create, admm_control_type, ddapc_control_create, ddapc_restraint_type, &
dft_control_create, dft_control_type, efield_type, qs_control_type, sccs_control_create, &
tddfpt2_control_type, tddfpt_control_create, tddfpt_control_type
tddfpt2_control_type, tddfpt_control_create, tddfpt_control_type, xtb_control_type
USE cp_files, ONLY: close_file,&
open_file
USE cp_log_handling, ONLY: cp_get_default_logger,&
Expand Down Expand Up @@ -780,7 +779,7 @@ SUBROUTINE read_qs_section(qs_control, qs_section)
CALL cite_reference(Seifert1996)
CASE (do_method_xtb)
qs_control%xtb = .TRUE.
!qxtb CALL cite_reference(Grimme2017)
CALL cite_reference(Grimme2017)
CASE (do_method_mndo)
CALL cite_reference(Dewar1977)
qs_control%semi_empirical = .TRUE.
Expand Down Expand Up @@ -1295,7 +1294,10 @@ SUBROUTINE write_dft_control(dft_control, dft_section)
IF (dft_control%qs_control%semi_empirical) RETURN
IF (dft_control%qs_control%dftb) RETURN
IF (dft_control%qs_control%xtb) RETURN
IF (dft_control%qs_control%xtb) THEN
CALL write_xtb_control(dft_control%qs_control%xtb_control, dft_section)
RETURN
END IF
CALL timeset(routineN, handle)
NULLIFY (logger)
Expand Down Expand Up @@ -1469,6 +1471,68 @@ SUBROUTINE write_dft_control(dft_control, dft_section)
END SUBROUTINE write_dft_control
! **************************************************************************************************
!> \brief Write the xTB control parameters to the output unit.
!> \param xtb_control ...
!> \param dft_section ...
! **************************************************************************************************
SUBROUTINE write_xtb_control(xtb_control, dft_section)
TYPE(xtb_control_type), POINTER :: xtb_control
TYPE(section_vals_type), POINTER :: dft_section
CHARACTER(len=*), PARAMETER :: routineN = 'write_xtb_control', &
routineP = moduleN//':'//routineN
INTEGER :: handle, output_unit
TYPE(cp_logger_type), POINTER :: logger
CALL timeset(routineN, handle)
NULLIFY (logger)
logger => cp_get_default_logger()
output_unit = cp_print_key_unit_nr(logger, dft_section, &
"PRINT%DFT_CONTROL_PARAMETERS", extension=".Log")
IF (output_unit > 0) THEN
WRITE (UNIT=output_unit, FMT="(/,T2,A,T31,A50)") &
"xTB| Parameter file", ADJUSTR(TRIM(xtb_control%parameter_file_name))
WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
"xTB| Basis expansion STO-NG", xtb_control%sto_ng
WRITE (UNIT=output_unit, FMT="(T2,A,T71,I10)") &
"xTB| Basis expansion STO-NG for Hydrogen", xtb_control%h_sto_ng
WRITE (UNIT=output_unit, FMT="(T2,A,T71,L10)") &
"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,T31,A50)") &
"xTB| D3 Dispersion: Parameter file", ADJUSTR(TRIM(xtb_control%dispersion_parameter_file))
WRITE (UNIT=output_unit, FMT="(T2,A,T51,3F10.3)") &
"xTB| Huckel constants ks kp kd", xtb_control%ks, xtb_control%kp, xtb_control%kd
WRITE (UNIT=output_unit, FMT="(T2,A,T61,2F10.3)") &
"xTB| Huckel constants ksp k2sh", xtb_control%ksp, xtb_control%k2sh
WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
"xTB| Mataga-Nishimoto exponent", xtb_control%kg
WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
"xTB| Repulsion potential exponent", xtb_control%kf
WRITE (UNIT=output_unit, FMT="(T2,A,T51,3F10.3)") &
"xTB| Coordination number scaling kcn(s) kcn(p) kcn(d)", &
xtb_control%kcns, xtb_control%kcnp, xtb_control%kcnd
WRITE (UNIT=output_unit, FMT="(T2,A,T71,F10.3)") &
"xTB| Electronegativity scaling", xtb_control%ken
WRITE (UNIT=output_unit, FMT="(T2,A,T61,2F10.3)") &
"xTB| Halogen potential scaling kxr kx2", xtb_control%kxr, xtb_control%kx2
WRITE (UNIT=output_unit, FMT="(/)")
END IF
CALL cp_print_key_finished_output(output_unit, logger, dft_section, &
"PRINT%DFT_CONTROL_PARAMETERS")
CALL timestop(handle)
END SUBROUTINE write_xtb_control
! **************************************************************************************************
!> \brief Purpose: Write the QS control parameters to the output unit.
!> \param qs_control ...
Expand Down
1 change: 1 addition & 0 deletions src/qs_environment.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -671,6 +671,7 @@ SUBROUTINE qs_init_subsys(qs_env, para_env, subsys, cell, cell_ref, use_ref_cell
zeff=qs_kind%xtb_parameter%zeff, zeff_correction=zeff_correction)
qs_kind%xtb_parameter%zeff = qs_kind%xtb_parameter%zeff-zeff_correction
END DO
!
! check for Ewald
IF (xtb_control%do_ewald) THEN
CALL ewald_env_create(ewald_env, para_env)
Expand Down
58 changes: 57 additions & 1 deletion src/qs_kind_types.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,9 @@ MODULE qs_kind_types
qs_control_type
USE cp_log_handling, ONLY: cp_get_default_logger,&
cp_logger_type
USE cp_output_handling, ONLY: cp_print_key_finished_output,&
USE cp_output_handling, ONLY: cp_p_file,&
cp_print_key_finished_output,&
cp_print_key_should_output,&
cp_print_key_unit_nr
USE cp_para_types, ONLY: cp_para_env_type
USE external_potential_types, ONLY: &
Expand Down Expand Up @@ -95,6 +97,7 @@ MODULE qs_kind_types
se_param_set_default
USE soft_basis_set, ONLY: create_soft_basis
USE string_utilities, ONLY: uppercase
USE xtb_parameters, ONLY: xtb_set_kab
USE xtb_types, ONLY: deallocate_xtb_atom_param,&
get_xtb_atom_param,&
write_xtb_atom_param,&
Expand Down Expand Up @@ -2483,12 +2486,65 @@ SUBROUTINE check_qs_kind_set(qs_kind_set, dft_control, subsys_section)
qs_kind => qs_kind_set(ikind)
CALL check_qs_kind(qs_kind, dft_control, subsys_section)
END DO
IF (dft_control%qs_control%xtb) THEN
CALL write_xtb_kab_param(qs_kind_set, subsys_section)
END IF
ELSE
CPABORT("The pointer qs_kind_set is not associated")
END IF
CALL timestop(handle)
END SUBROUTINE check_qs_kind_set

! **************************************************************************************************
!> \brief ...
!> \param qs_kind_set ...
!> \param subsys_section ...
! **************************************************************************************************
SUBROUTINE write_xtb_kab_param(qs_kind_set, subsys_section)

TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
TYPE(section_vals_type), POINTER :: subsys_section

CHARACTER(len=*), PARAMETER :: routineN = 'write_xtb_kab_param', &
routineP = moduleN//':'//routineN

CHARACTER(LEN=default_string_length) :: aname, bname
INTEGER :: ikind, io_unit, jkind, nkind, za, zb
TYPE(cp_logger_type), POINTER :: logger
TYPE(qs_kind_type), POINTER :: qs_kinda, qs_kindb
TYPE(xtb_atom_type), POINTER :: xtb_parameter_a, xtb_parameter_b

NULLIFY (logger)
logger => cp_get_default_logger()
IF (BTEST(cp_print_key_should_output(logger%iter_info, subsys_section, &
"PRINT%KINDS/POTENTIAL"), cp_p_file)) THEN

io_unit = cp_print_key_unit_nr(logger, subsys_section, "PRINT%KINDS", extension=".Log")
IF (io_unit > 0) THEN

WRITE (io_unit, "(/,T2,A)") "xTB| Kab parameters"
nkind = SIZE(qs_kind_set)
DO ikind = 1, nkind
qs_kinda => qs_kind_set(ikind)
CALL get_qs_kind(qs_kinda, xtb_parameter=xtb_parameter_a)
CALL get_xtb_atom_param(xtb_parameter_a, aname=aname, z=za)
DO jkind = ikind, nkind
qs_kindb => qs_kind_set(jkind)
CALL get_qs_kind(qs_kindb, xtb_parameter=xtb_parameter_b)
CALL get_xtb_atom_param(xtb_parameter_b, aname=bname, z=zb)
WRITE (io_unit, "(A,T10,A15,T25,A15,T71,F10.3)") &
" Kab:", TRIM(aname), TRIM(bname), xtb_set_kab(za, zb)
END DO
END DO
WRITE (io_unit, *)

END IF

CALL cp_print_key_finished_output(io_unit, logger, subsys_section, "PRINT%KINDS")
END IF

END SUBROUTINE write_xtb_kab_param

! **************************************************************************************************
!> \brief Set the components of an atomic kind data set.
!> \param qs_kind ...
Expand Down
29 changes: 23 additions & 6 deletions src/xtb_types.F
View file
Original file line number Diff line number Diff line change
Expand Up @@ -311,10 +311,13 @@ SUBROUTINE write_xtb_atom_param(xtb_parameter, subsys_section)
CHARACTER(LEN=*), PARAMETER :: routineN = 'write_xtb_atom_param', &
routineP = moduleN//':'//routineN

CHARACTER(LEN=default_string_length) :: aname, typ
INTEGER :: io_unit, lmax, natorb, z
CHARACTER(LEN=default_string_length) :: aname, bb
INTEGER :: i, io_unit, m, natorb, nshell
INTEGER, DIMENSION(5) :: lval, nval, occupation
LOGICAL :: defined
REAL(dp) :: zeff
REAL(KIND=dp) :: alpha, en, eta, xgamma, zneff
REAL(KIND=dp), DIMENSION(5) :: hen, kappa, kpoly, zeta
TYPE(cp_logger_type), POINTER :: logger

NULLIFY (logger)
Expand All @@ -327,19 +330,33 @@ SUBROUTINE write_xtb_atom_param(xtb_parameter, subsys_section)
extension=".Log")

IF (io_unit > 0) THEN
CALL get_xtb_atom_param(xtb_parameter, aname=aname, typ=typ, defined=defined, &
z=z, zeff=zeff, natorb=natorb, lmax=lmax)
CALL get_xtb_atom_param(xtb_parameter, aname=aname, defined=defined, zeff=zeff, natorb=natorb)
CALL get_xtb_atom_param(xtb_parameter, nshell=nshell, lval=lval, nval=nval, occupation=occupation)
CALL get_xtb_atom_param(xtb_parameter, kpoly=kpoly, kappa=kappa, hen=hen, zeta=zeta)
CALL get_xtb_atom_param(xtb_parameter, electronegativity=en, xgamma=xgamma, eta=eta, alpha=alpha, zneff=zneff)

bb = " "
WRITE (UNIT=io_unit, FMT="(/,A,T67,A14)") " xTB parameters: ", TRIM(aname)
IF (defined) THEN
m = 5-nshell
WRITE (UNIT=io_unit, FMT="(T16,A,T71,F10.2)") "Effective core charge:", zeff
WRITE (UNIT=io_unit, FMT="(T16,A,T71,I10)") "Number of orbitals:", natorb
WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5(A4,I1,I2,A1))") "Basis set [nl]", bb(1:8*m), &
(" [", nval(i), lval(i), "]", i=1, nshell)
WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5F8.3)") "Slater Exponent", bb(1:8*m), (zeta(i), i=1, nshell)
WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5I8)") "Ref. occupation", bb(1:8*m), (occupation(i), i=1, nshell)
WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5F8.3)") "Energy levels [au]", bb(1:8*m), (hen(i), i=1, nshell)
WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5F8.3)") "Kpoly", bb(1:8*m), (kpoly(i), i=1, nshell)
WRITE (UNIT=io_unit, FMT="(T16,A,T71,F10.3)") "Electronegativity", en
WRITE (UNIT=io_unit, FMT="(T16,A,T71,F10.3)") "Mataga-Nishimoto constant (eta)", eta
WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5F8.3)") "Mataga-Nishimoto scaling kappa", bb(1:8*m), (kappa(i), i=1, nshell)
WRITE (UNIT=io_unit, FMT="(T16,A,T71,F10.3)") "3rd Order constant", xgamma
WRITE (UNIT=io_unit, FMT="(T16,A,T61,2F10.3)") "Repulsion potential [Z,alpha]", zneff, alpha
ELSE
WRITE (UNIT=io_unit, FMT="(T55,A)") "Parameters are not defined"
END IF
END IF
CALL cp_print_key_finished_output(io_unit, logger, subsys_section, &
"PRINT%KINDS")
CALL cp_print_key_finished_output(io_unit, logger, subsys_section, "PRINT%KINDS")
END IF

END SUBROUTINE write_xtb_atom_param
Expand Down

0 comments on commit e5cf0c8

Please sign in to comment.