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

Quickstep: STO-NG basis sets

Browse source 
svn-origin-rev: 18583
This commit is contained in:
Jürg Hutter 2018-08-10 16:07:01 +00:00
parent bebe6d1c4f
commit 295087e08b
8 changed files with 331 additions and 11 deletions
Download patch file Download diff file Expand all files Collapse all files

174
src/aobasis/basis_set_types.F
View file

@ -123,6 +123,7 @@ MODULE basis_set_types
write_orb_basis_set
PUBLIC :: allocate_sto_basis_set, &
read_sto_basis_set, &
create_gto_from_sto_basis, &
deallocate_sto_basis_set, &
set_sto_basis_set, &
@ -2223,6 +2224,179 @@ CONTAINS
END SUBROUTINE set_sto_basis_set
! **************************************************************************************************
!> \brief ...
!> \param element_symbol ...
!> \param basis_set_name ...
!> \param sto_basis_set ...
!> \param para_env ...
!> \param dft_section ...
! **************************************************************************************************
SUBROUTINE read_sto_basis_set(element_symbol, basis_set_name, sto_basis_set, para_env, dft_section)
! Read a Slater-type orbital (STO) basis set from the database file.
CHARACTER(LEN=*), INTENT(IN) :: element_symbol, basis_set_name
TYPE(sto_basis_set_type), POINTER :: sto_basis_set
TYPE(cp_para_env_type), POINTER :: para_env
TYPE(section_vals_type), POINTER :: dft_section
CHARACTER(len=*), PARAMETER :: routineN = 'read_sto_basis_set', &
routineP = moduleN//':'//routineN
CHARACTER(LEN=10) :: nlsym
CHARACTER(LEN=2) :: lsym
CHARACTER(LEN=240) :: line
CHARACTER(LEN=242) :: line2
CHARACTER(len=default_path_length) :: basis_set_file_name, tmp
CHARACTER(LEN=default_path_length), DIMENSION(:), &
POINTER :: cbasis
CHARACTER(LEN=LEN(basis_set_name)) :: bsname
CHARACTER(LEN=LEN(basis_set_name)+2) :: bsname2
CHARACTER(LEN=LEN(element_symbol)) :: symbol
CHARACTER(LEN=LEN(element_symbol)+2) :: symbol2
INTEGER :: ibasis, irep, iset, nbasis, nq, nset, &
strlen1, strlen2
LOGICAL :: basis_found, found, match
REAL(KIND=dp) :: zet
TYPE(cp_parser_type), POINTER :: parser
line = ""
line2 = ""
symbol = ""
symbol2 = ""
bsname = ""
bsname2 = ""
nbasis = 1
sto_basis_set%name = basis_set_name
CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
n_rep_val=nbasis)
ALLOCATE (cbasis(nbasis))
DO ibasis = 1, nbasis
CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
i_rep_val=ibasis, c_val=cbasis(ibasis))
NULLIFY (parser)
basis_set_file_name = cbasis(ibasis)
tmp = basis_set_file_name
CALL uppercase(tmp)
END DO
! Search for the requested basis set in the basis set file
! until the basis set is found or the end of file is reached
basis_found = .FALSE.
basis_loop: DO ibasis = 1, nbasis
IF (basis_found) EXIT basis_loop
NULLIFY (parser)
basis_set_file_name = cbasis(ibasis)
CALL parser_create(parser, basis_set_file_name, para_env=para_env)
bsname = basis_set_name
symbol = element_symbol
irep = 0
tmp = basis_set_name
CALL uppercase(tmp)
IF (tmp .NE. "NONE") THEN
search_loop: DO
CALL parser_search_string(parser, TRIM(bsname), .TRUE., found, line)
IF (found) THEN
CALL uppercase(symbol)
CALL uppercase(bsname)
match = .FALSE.
CALL uppercase(line)
! Check both the element symbol and the basis set name
line2 = " "//line//" "
symbol2 = " "//TRIM(symbol)//" "
bsname2 = " "//TRIM(bsname)//" "
strlen1 = LEN_TRIM(symbol2)+1
strlen2 = LEN_TRIM(bsname2)+1
IF ((INDEX(line2, symbol2(:strlen1)) > 0) .AND. &
(INDEX(line2, bsname2(:strlen2)) > 0)) match = .TRUE.
IF (match) THEN
! Read the basis set information
CALL parser_get_object(parser, nset, newline=.TRUE.)
sto_basis_set%nshell = nset
CALL reallocate(sto_basis_set%nq, 1, nset)
CALL reallocate(sto_basis_set%lq, 1, nset)
CALL reallocate(sto_basis_set%zet, 1, nset)
ALLOCATE (sto_basis_set%symbol(nset))
DO iset = 1, nset
CALL parser_get_object(parser, nq, newline=.TRUE.)
CALL parser_get_object(parser, lsym)
CALL parser_get_object(parser, zet)
sto_basis_set%nq(iset) = nq
sto_basis_set%zet(iset) = zet
WRITE (nlsym, "(I2,A)") nq, TRIM(lsym)
sto_basis_set%symbol(iset) = TRIM(nlsym)
SELECT CASE (TRIM (lsym))
CASE ("S", "s")
sto_basis_set%lq(iset) = 0
CASE ("P", "p")
sto_basis_set%lq(iset) = 1
CASE ("D", "d")
sto_basis_set%lq(iset) = 2
CASE ("F", "f")
sto_basis_set%lq(iset) = 3
CASE ("G", "g")
sto_basis_set%lq(iset) = 4
CASE ("H", "h")
sto_basis_set%lq(iset) = 5
CASE ("I", "i", "J", "j")
sto_basis_set%lq(iset) = 6
CASE ("K", "k")
sto_basis_set%lq(iset) = 7
CASE ("L", "l")
sto_basis_set%lq(iset) = 8
CASE ("M", "m")
sto_basis_set%lq(iset) = 9
CASE DEFAULT
CALL cp_abort(__LOCATION__, &
"The requested basis set <"//TRIM(bsname)// &
"> for element <"//TRIM(symbol)//"> has an ivalid component: ")
END SELECT
END DO
basis_found = .TRUE.
EXIT search_loop
END IF
ELSE
EXIT search_loop
END IF
END DO search_loop
ELSE
match = .FALSE.
ENDIF
CALL parser_release(parser)
END DO basis_loop
IF (tmp .NE. "NONE") THEN
IF (.NOT. basis_found) THEN
basis_set_file_name = ""
DO ibasis = 1, nbasis
basis_set_file_name = TRIM(basis_set_file_name)//"<"//TRIM(cbasis(ibasis))//"> "
END DO
CALL cp_abort(__LOCATION__, &
"The requested basis set <"//TRIM(bsname)// &
"> for element <"//TRIM(symbol)//"> was not "// &
"found in the basis set files "// &
TRIM(basis_set_file_name))
END IF
END IF
DEALLOCATE (cbasis)
END SUBROUTINE read_sto_basis_set
! **************************************************************************************************
!> \brief ...
!> \param sto_basis_set ...

6
src/input_cp2k_dft.F
View file

@ -2961,6 +2961,12 @@ CONTAINS
CALL keyword_release(keyword)
! Integers
CALL keyword_create(keyword, name="STO_NG", &
description="Order of Gaussian type expansion of Slater orbital basis sets.", &
usage="STO_NG", default_i_val=6)
CALL section_add_keyword(section, keyword)
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="LMAXN1", &
variants=(/"LMAXRHO1"/), &
description="GAPW : max L number for expansion of the atomic densities in spherical gaussians", &

7
src/input_cp2k_subsys.F
View file

@ -1064,8 +1064,11 @@ CONTAINS
CALL keyword_release(keyword)
CALL keyword_create(keyword, name="BASIS_SET", &
description="The primary Gaussian basis set (NONE implies no basis used, meaningful with GHOST)", &
usage="BASIS_SET [type] DZVP", type_of_var=char_t, default_c_vals=(/" ", " "/), &
description="The primary Gaussian basis set (NONE implies no basis used, meaningful with GHOST)"// &
"Defaults are set for TYPE {ORB} and FORM {GTO}. Possible values for TYPE are {ORB, AUX, RI_AUX, LRI, ...}. "// &
"Possible values for FORM are {GTO, STO}. Where STO results in a GTO expansion of a Slater type basis."// &
"If a value for FORM is given, also TYPE has to be set explicitely.", &
usage="BASIS_SET [type] [form] DZVP", type_of_var=char_t, default_c_vals=(/" ", " ", " "/), &
citations=(/VandeVondele2005a, VandeVondele2007/), &
repeats=.TRUE., n_var=-1)
CALL section_add_keyword(section, keyword)

52
src/qs_kind_types.F
View file

@ -35,9 +35,11 @@ MODULE qs_kind_types
remove_basis_from_container,&
remove_basis_set_container
USE basis_set_types, ONLY: &
allocate_gto_basis_set, combine_basis_sets, deallocate_gto_basis_set, get_gto_basis_set, &
gto_basis_set_type, init_aux_basis_set, init_orb_basis_set, read_gto_basis_set, &
write_gto_basis_set, write_orb_basis_set
allocate_gto_basis_set, allocate_sto_basis_set, combine_basis_sets, &
create_gto_from_sto_basis, deallocate_gto_basis_set, deallocate_sto_basis_set, &
get_gto_basis_set, gto_basis_set_type, init_aux_basis_set, init_orb_basis_set, &
read_gto_basis_set, read_sto_basis_set, sto_basis_set_type, write_gto_basis_set, &
write_orb_basis_set
USE cp_control_types, ONLY: dft_control_type,&
qs_control_type
USE cp_log_handling, ONLY: cp_get_default_logger,&
@ -1372,7 +1374,8 @@ CONTAINS
CHARACTER(LEN=default_string_length), &
DIMENSION(:), POINTER :: tmpstringlist
CHARACTER(LEN=default_string_length), &
DIMENSION(maxbas) :: basis_set_name, basis_set_type
DIMENSION(maxbas) :: basis_set_form, basis_set_name, &
basis_set_type
INTEGER :: handle, i, i_rep, ipaodesc, ipaopot, ipos, j, jj, k_rep, l, m, n_rep, nb_rep, &
nexp, ngauss, nlcc, nloc, nnl, norbitals, npaodesc, npaopot, nppnl, nspin, z
INTEGER, DIMENSION(:), POINTER :: add_el, elec_conf, orbitals
@ -1390,6 +1393,7 @@ CONTAINS
TYPE(section_vals_type), POINTER :: basis_section, bs_section, dft_plus_u_section, &
dft_section, enforce_occupation_section, kgpot_section, pao_desc_section, &
pao_pot_section, potential_section, spin_section
TYPE(sto_basis_set_type), POINTER :: sto_basis_set
CALL timeset(routineN, handle)
@ -1398,6 +1402,7 @@ CONTAINS
update_input = .TRUE.
basis_set_name(:) = ""
basis_set_type(:) = ""
basis_set_form(:) = ""
aux_basis_set_name = ""
ri_aux_basis_set_name = ""
orb_basis_set_name = ""
@ -1494,15 +1499,26 @@ CONTAINS
CALL section_vals_val_get(kind_section, i_rep_section=k_rep, &
keyword_name="BASIS_SET", i_rep_val=i, c_vals=tmpstringlist)
IF (SIZE(tmpstringlist) == 1) THEN
! default is orbital type
! default is orbital type and GTO
basis_set_type(i) = "ORB"
basis_set_form(i) = "GTO"
basis_set_name(i) = tmpstringlist(1)
ELSEIF (SIZE(tmpstringlist) == 2) THEN
! default is GTO
basis_set_type(i) = tmpstringlist(1)
basis_set_form(i) = "GTO"
basis_set_name(i) = tmpstringlist(2)
ELSEIF (SIZE(tmpstringlist) == 3) THEN
basis_set_type(i) = tmpstringlist(1)
basis_set_form(i) = tmpstringlist(2)
basis_set_name(i) = tmpstringlist(3)
ELSE
CPABORT("")
END IF
! check that we have a valid basis set form
IF (basis_set_form(i) /= "GTO" .AND. basis_set_form(i) /= "STO") THEN
CPABORT("BASIS_SET_FORM invalid")
END IF
END DO
! Deprecated keywords
@ -1540,6 +1556,7 @@ CONTAINS
nb_rep = nb_rep+1
CPASSERT(nb_rep <= maxbas)
basis_set_type(nb_rep) = "RI_AUX"
basis_set_form(nb_rep) = "GTO"
basis_set_name(nb_rep) = ri_aux_basis_set_name
END IF
! end of deprecated input
@ -1975,13 +1992,30 @@ CONTAINS
CASE DEFAULT
! set nauss for STO expansion
CALL section_vals_val_get(dft_section, "QS%STO_NG", i_val=ngauss)
! Allocate and initialise the basis set data set structure
! first external basis sets
DO i = 1, nb_rep
NULLIFY (tmp_basis_set)
CALL allocate_gto_basis_set(tmp_basis_set)
CALL read_gto_basis_set(qs_kind%element_symbol, basis_set_name(i), &
tmp_basis_set, para_env, dft_section)
SELECT CASE (basis_set_form (i))
CASE ("GTO")
NULLIFY (tmp_basis_set)
CALL allocate_gto_basis_set(tmp_basis_set)
CALL read_gto_basis_set(qs_kind%element_symbol, basis_set_name(i), &
tmp_basis_set, para_env, dft_section)
CASE ("STO")
NULLIFY (sto_basis_set)
CALL allocate_sto_basis_set(sto_basis_set)
CALL read_sto_basis_set(qs_kind%element_symbol, basis_set_name(i), &
sto_basis_set, para_env, dft_section)
NULLIFY (tmp_basis_set)
CALL create_gto_from_sto_basis(sto_basis_set, tmp_basis_set, ngauss)
CALL deallocate_sto_basis_set(sto_basis_set)
CASE DEFAULT
CALL cp_abort(__LOCATION__, &
"Invalid basis set form "//TRIM(basis_set_form(i))// &
"for atomic kind <"//TRIM(qs_kind%name)//">")
END SELECT
tmp = basis_set_type(i)
CALL uppercase(tmp)
CALL add_basis_set_to_container(qs_kind%basis_sets, tmp_basis_set, tmp)

44
tests/QS/regtest-sto/BASIS_STO Normal file
View file

@ -0,0 +1,44 @@
H DZ_AE
2
1 S 0.760
1 S 1.280
H QZ_AE
8
1 S 3.300
1 S 2.000
1 S 1.400
1 S 1.000
1 S 0.710
2 P 2.000
2 P 1.000
3 D 2.500
3 D 1.500
O DZ_AE
6
1 S 9.80
1 S 6.80
2 S 1.70
2 S 2.82
2 P 3.06
2 P 1.30
O QZ_AE
16
1 S 14.800
1 S 8.750
1 S 6.650
2 S 7.100
2 S 3.200
2 S 2.050
2 S 1.500
2 S 0.750
2 P 5.700
2 P 3.050
2 P 1.650
2 P 1.000
3 D 2.500
3 D 1.500
4 F 4.000
4 F 2.000

50
tests/QS/regtest-sto/H2O_t1.inp Normal file
View file

@ -0,0 +1,50 @@
&FORCE_EVAL
METHOD Quickstep
&DFT
BASIS_SET_FILE_NAME ./BASIS_STO
POTENTIAL_FILE_NAME POTENTIAL
&MGRID
CUTOFF 250
NGRIDS 5
&END MGRID
&QS
METHOD GAPW
STO_NG 4
&END QS
&SCF
EPS_SCF 1.0E-6
SCF_GUESS ATOMIC
&OT ON
MINIMIZER DIIS
PRECONDITIONER FULL_ALL
&END OT
&END SCF
&XC
&XC_FUNCTIONAL Pade
&END XC_FUNCTIONAL
&END XC
&END DFT
&SUBSYS
&CELL
ABC 6.0 6.0 6.0
&END CELL
&COORD
O 0.000000 0.000000 -0.065587
H 0.000000 -0.757136 0.520545
H 0.000000 0.757136 0.520545
&END COORD
&KIND H
BASIS_SET ORB STO DZ_AE
POTENTIAL ALL
&END KIND
&KIND O
BASIS_SET ORB STO DZ_AE
POTENTIAL ALL
&END KIND
&END SUBSYS
&END FORCE_EVAL
&GLOBAL
PROJECT sto1
PRINT_LEVEL MEDIUM
&END GLOBAL

8
tests/QS/regtest-sto/TEST_FILES Normal file
View file

@ -0,0 +1,8 @@
# runs are executed in the same order as in this file
# the second field tells which test should be run in order to compare with the last available output
# e.g. 0 means do not compare anything, running is enough
# 1 compares the last total energy in the file
# for details see cp2k/tools/do_regtest
# tests GAPW STO
H2O_t1.inp 1 1e-13 -75.659689099442119
#EOF

1
tests/QS/regtest-sto/TEST_FILES_RESET Normal file
View file

@ -0,0 +1 @@
#