Skip to content

Commit

Permalink
Remove auto basis option for ADMM (cp2k#3307)
Browse files Browse the repository at this point in the history
  • Loading branch information
juerghutter authored Mar 8, 2024
1 parent 2204749 commit 8976884
Show file tree
Hide file tree
Showing 3 changed files with 16 additions and 227 deletions.
218 changes: 1 addition & 217 deletions src/auto_basis.F
Original file line number Diff line number Diff line change
Expand Up @@ -28,8 +28,6 @@ MODULE auto_basis
rootpi
USE orbital_pointers, ONLY: init_orbital_pointers
USE periodic_table, ONLY: get_ptable_info
USE powell, ONLY: opt_state_type,&
powell_optimize
USE qs_kind_types, ONLY: get_qs_kind,&
qs_kind_type
#include "./base/base_uses.f90"
Expand All @@ -40,7 +38,7 @@ MODULE auto_basis

CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'auto_basis'

PUBLIC :: create_ri_aux_basis_set, create_aux_fit_basis_set, create_lri_aux_basis_set, &
PUBLIC :: create_ri_aux_basis_set, create_lri_aux_basis_set, &
create_oce_basis

CONTAINS
Expand Down Expand Up @@ -395,220 +393,6 @@ SUBROUTINE create_lri_aux_basis_set(lri_aux_basis_set, qs_kind, basis_cntrl, &
END IF

END SUBROUTINE create_lri_aux_basis_set
! **************************************************************************************************
!> \brief Create a AUX_FIT basis set using some heuristics
!> \param aux_fit_basis ...
!> \param qs_kind ...
!> \param basis_cntrl ...
!> \date 01.11.2017
!> \author JGH
! **************************************************************************************************
SUBROUTINE create_aux_fit_basis_set(aux_fit_basis, qs_kind, basis_cntrl)
TYPE(gto_basis_set_type), POINTER :: aux_fit_basis
TYPE(qs_kind_type), INTENT(IN) :: qs_kind
INTEGER, INTENT(IN) :: basis_cntrl

CHARACTER(LEN=2) :: element_symbol
CHARACTER(LEN=default_string_length) :: bsname
INTEGER :: i, iset, l, laux, lmax, lval, maxpgf, &
mx, nf, np, nsets, nx, z
INTEGER, DIMENSION(0:9) :: nfun, nval
INTEGER, DIMENSION(0:9, 1:20) :: nl
INTEGER, DIMENSION(1:20) :: lset, npgf
INTEGER, DIMENSION(:), POINTER :: econf
REAL(KIND=dp) :: amet, z1, z2, zvel
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: zval
REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: gcval, zet
REAL(KIND=dp), DIMENSION(0:9) :: zeff, zmax, zmin
REAL(KIND=dp), DIMENSION(25) :: afit, xval
TYPE(gto_basis_set_type), POINTER :: orb_basis_set
TYPE(opt_state_type) :: ostate

!
CPABORT("Automatic basis set generation not activated")
!
CPASSERT(.NOT. ASSOCIATED(aux_fit_basis))
NULLIFY (orb_basis_set)
CALL get_qs_kind(qs_kind, basis_set=orb_basis_set, basis_type="ORB")
IF (ASSOCIATED(orb_basis_set)) THEN
CALL get_basis_keyfigures(orb_basis_set, lmax, zmin, zmax, zeff)
CALL get_qs_kind(qs_kind, zeff=zvel, elec_conf=econf, element_symbol=element_symbol)
CALL get_ptable_info(element_symbol, ielement=z)
lval = 0
DO l = 0, MAXVAL(UBOUND(econf))
IF (econf(l) > 0) lval = l
END DO
IF (SUM(econf) /= NINT(zvel)) THEN
CPWARN("Valence charge and electron configuration not consistent")
END IF
nval = 0
DO l = 0, lval
nx = econf(l)
DO
IF (nx > 0) THEN
nval(l) = nval(l) + 1
nx = nx - 2*(2*l + 1)
ELSE
EXIT
END IF
END DO
END DO
nfun = nval
SELECT CASE (basis_cntrl)
CASE (0)
laux = lval
DO l = 0, lval
nfun(l) = nfun(l) + 1
END DO
CASE (1)
laux = MIN(lval + 1, lmax)
DO l = 0, lval
nfun(l) = nfun(l) + 1
END DO
IF (laux > lval) nfun(laux) = 1
CASE (2)
laux = MIN(lval + 1, lmax)
DO l = 0, lval
nfun(l) = nfun(l) + 2
END DO
IF (laux > lval) nfun(laux) = 1
CASE (3)
laux = MIN(lval + 2, lmax)
DO l = 0, lval
nfun(l) = nfun(l) + 3
END DO
IF (laux > lval) nfun(lval + 1) = 2
IF (laux > lval + 1) nfun(laux) = 1
CASE DEFAULT
CPABORT("Invalid value of control variable")
END SELECT
!
nsets = 0
maxpgf = 0
DO l = 0, lval
z1 = MAX(zmin(l), 0.10_dp + l*0.025_dp)
z2 = zmax(l)
mx = CEILING(LOG(z2/z1)/LOG(5.0_dp))
IF (nval(l) > 1) THEN
nsets = nsets + 2
maxpgf = MAX(maxpgf, nval(l), mx, 3)
ELSEIF (nval(l) == 1) THEN
nsets = nsets + 1
maxpgf = MAX(maxpgf, mx, 3)
END IF
DO i = nval(l) + 1, nfun(l)
maxpgf = MAX(maxpgf, mx, 1)
nsets = nsets + 1
END DO
END DO
DO l = lval + 1, laux
maxpgf = MAX(maxpgf, 1)
nsets = nsets + nfun(l)
END DO
!
ALLOCATE (zet(maxpgf, nsets))
zet = 0.0_dp
nl = 0
iset = 0
DO l = 0, laux
amet = 2.50_dp
! optimize exponensts
z1 = MAX(zmin(l), 0.20_dp + l*0.025_dp)
z2 = zmax(l)
mx = CEILING(LOG(z2/z1)/LOG(4.0_dp))
IF (nval(l) > 1) THEN
nx = MAX(nfun(l), mx, 3)
ELSE IF (nval(l) == 1) THEN
nx = MAX(nfun(l), mx, 2)
ELSE
nx = nfun(l)
END IF
CPASSERT(nx <= 25)
! Powell optimizer
CALL get_basis_functions(orb_basis_set, l, np, nf, zval, gcval)
ostate%nf = 0
ostate%nvar = nx
xval(:) = 5.00_dp
xval(1) = z1
ostate%rhoend = 1.e-8_dp
ostate%rhobeg = 5.e-1_dp
ostate%maxfun = 1000
ostate%iprint = 0
ostate%unit = 0
ostate%state = 0
DO
IF (ostate%state == 2) THEN
afit(1) = xval(1)
DO i = 2, nx
afit(i) = afit(i - 1)*xval(i)
END DO
CALL overlap_maximum(l, np, nf, zval, gcval, nx, afit, amet, ostate%f)
CALL neb_potential(xval, ostate%nvar, ostate%f)
END IF
IF (ostate%state == -1) EXIT
CALL powell_optimize(ostate%nvar, xval, ostate)
END DO
ostate%state = 8
CALL powell_optimize(ostate%nvar, xval, ostate)
afit(1) = xval(1)
DO i = 2, nx
afit(i) = afit(i - 1)*xval(i)
END DO
DEALLOCATE (zval, gcval)
!
IF (nval(l) > 1) THEN
! split set
iset = iset + 1
lset(iset) = l
npgf(iset) = nx - 1
nl(l, iset) = nval(l) - 1
zet(1:nx - 1, iset) = afit(1:nx - 1)
! new set
iset = iset + 1
lset(iset) = l
npgf(iset) = nx - 1
nl(l, iset) = 1
zet(1:nx - 1, iset) = afit(2:nx)
DO i = 1, nfun(l) - 2
iset = iset + 1
lset(iset) = l
npgf(iset) = 1
zet(1, iset) = afit(nx - i + 1)
nl(l, iset) = 1
END DO
ELSEIF (nval(l) == 1) THEN
iset = iset + 1
lset(iset) = l
npgf(iset) = nx
zet(1:nx, iset) = afit(1:nx)
nl(l, iset) = 1
DO i = 1, nfun(l) - 1
iset = iset + 1
lset(iset) = l
npgf(iset) = 1
zet(1, iset) = afit(nx - i + 1)
nl(l, iset) = 1
END DO
ELSE
DO i = 1, nfun(l)
iset = iset + 1
lset(iset) = l
npgf(iset) = 1
zet(1, iset) = afit(i)
nl(l, iset) = 1
END DO
END IF
END DO
!
bsname = TRIM(element_symbol)//"-AUX-FIT-"//TRIM(orb_basis_set%name)
!
CALL create_aux_basis(aux_fit_basis, bsname, nsets, lset, lset, nl, npgf, zet)
!
DEALLOCATE (zet)
!
END IF

END SUBROUTINE create_aux_fit_basis_set

! **************************************************************************************************
!> \brief ...
Expand Down
12 changes: 11 additions & 1 deletion src/min_basis_set.F
Original file line number Diff line number Diff line change
Expand Up @@ -140,7 +140,7 @@ SUBROUTINE create_min_basis(min_basis_set, zval, econf, mao, ngau)
CHARACTER(len=6) :: str
CHARACTER(len=6), DIMENSION(:), POINTER :: sym
CHARACTER(LEN=default_string_length) :: bname
INTEGER :: i, iss, l, n, nmao, nn, nss
INTEGER :: i, iss, l, lm, n, nmao, nn, nss
INTEGER, DIMENSION(0:3) :: nae, nco, npe
INTEGER, DIMENSION(4, 7) :: ne
INTEGER, DIMENSION(:), POINTER :: lq, nq
Expand Down Expand Up @@ -186,6 +186,7 @@ SUBROUTINE create_min_basis(min_basis_set, zval, econf, mao, ngau)
CASE (2)
npe(0) = npe(0) + 2
CASE (3)
npe(1) = npe(1) + 1
CASE (4)
npe(0) = npe(0) + 1
npe(1) = npe(1) + 1
Expand Down Expand Up @@ -215,6 +216,15 @@ SUBROUTINE create_min_basis(min_basis_set, zval, econf, mao, ngau)
CALL cp_warn(__LOCATION__, "Reference basis has been adjusted according to MAO value.")
END IF

! All shells should have at least 1 function
lm = 0
DO l = 0, 3
IF (npe(l) > 0) lm = l
END DO
DO l = 0, lm
IF (npe(l) == 0) npe(l) = 1
END DO

nss = SUM(npe)
ALLOCATE (sym(nss), lq(nss), nq(nss), zet(nss))
iss = 0
Expand Down
13 changes: 4 additions & 9 deletions src/qs_environment.F
Original file line number Diff line number Diff line change
Expand Up @@ -20,8 +20,7 @@ MODULE qs_environment
USE almo_scf_env_methods, ONLY: almo_scf_env_create
USE atom_kind_orbitals, ONLY: calculate_atomic_relkin
USE atomic_kind_types, ONLY: atomic_kind_type
USE auto_basis, ONLY: create_aux_fit_basis_set,&
create_lri_aux_basis_set,&
USE auto_basis, ONLY: create_lri_aux_basis_set,&
create_ri_aux_basis_set
USE basis_set_container_types, ONLY: add_basis_set_to_container
USE basis_set_types, ONLY: basis_sort_zet,&
Expand Down Expand Up @@ -840,19 +839,15 @@ SUBROUTINE qs_init_subsys(qs_env, para_env, subsys, cell, cell_ref, use_ref_cell
CALL get_qs_kind_set(qs_kind_set, dft_plus_u_atom_present=dft_control%dft_plus_u)

IF (dft_control%do_admm) THEN
! Check if ADMM basis is available, auto-generate if needed
! Check if ADMM basis is available
CALL get_qs_env(qs_env, nkind=nkind)
DO ikind = 1, nkind
NULLIFY (aux_fit_basis)
qs_kind => qs_kind_set(ikind)
CALL get_qs_kind(qs_kind, basis_set=aux_fit_basis, basis_type="AUX_FIT")
IF (.NOT. (ASSOCIATED(aux_fit_basis))) THEN
! AUX_FIT basis set is not yet loaded
CALL cp_warn(__LOCATION__, "Automatic Generation of AUX_FIT basis. "// &
"This is experimental code.")
! Generate a default basis
CALL create_aux_fit_basis_set(aux_fit_basis, qs_kind, dft_control%auto_basis_aux_fit)
CALL add_basis_set_to_container(qs_kind%basis_sets, aux_fit_basis, "AUX_FIT")
! AUX_FIT basis set is not available
CPABORT("AUX_FIT basis set is not defined. ")
END IF
END DO
END IF
Expand Down

0 comments on commit 8976884

Please sign in to comment.