1!--------------------------------------------------------------------------------------------------! 2! CP2K: A general program to perform molecular dynamics simulations ! 3! Copyright (C) 2000 - 2020 CP2K developers group ! 4!--------------------------------------------------------------------------------------------------! 5 6! ************************************************************************************************** 7!> \brief Print basis sets in CP2K format 8!> \par History 9!> \author JGH (12.2017) 10! ************************************************************************************************** 11MODULE basis_set_output 12 USE basis_set_types, ONLY: get_gto_basis_set,& 13 gto_basis_set_type 14 USE cp2k_info, ONLY: compile_revision,& 15 cp2k_version,& 16 r_datx,& 17 r_host_name,& 18 r_user_name 19 USE cp_files, ONLY: close_file,& 20 open_file 21 USE cp_log_handling, ONLY: cp_get_default_logger,& 22 cp_logger_get_default_io_unit,& 23 cp_logger_type 24 USE input_section_types, ONLY: section_vals_type,& 25 section_vals_val_get 26 USE kinds, ONLY: default_string_length,& 27 dp 28 USE qs_environment_types, ONLY: get_qs_env,& 29 qs_environment_type 30 USE qs_kind_types, ONLY: get_qs_kind,& 31 qs_kind_type 32#include "./base/base_uses.f90" 33 34 IMPLICIT NONE 35 PRIVATE 36 37 ! Global parameters 38 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'basis_set_output' 39 PUBLIC :: print_basis_set_file 40 41! ************************************************************************************************** 42 43CONTAINS 44 45! ************************************************************************************************** 46!> \brief ... 47!> \param qs_env ... 48!> \param base_section ... 49! ************************************************************************************************** 50 SUBROUTINE print_basis_set_file(qs_env, base_section) 51 52 TYPE(qs_environment_type), POINTER :: qs_env 53 TYPE(section_vals_type), POINTER :: base_section 54 55 CHARACTER(len=*), PARAMETER :: routineN = 'print_basis_set_file', & 56 routineP = moduleN//':'//routineN 57 58 CHARACTER(LEN=2) :: element_symbol 59 CHARACTER(LEN=default_string_length) :: bname, filename 60 INTEGER :: ikind, iunit, nkind, ounit 61 INTEGER, SAVE :: ncalls = 0 62 TYPE(cp_logger_type), POINTER :: logger 63 TYPE(gto_basis_set_type), POINTER :: aux_fit_basis, lri_aux_basis, orb_basis, & 64 ri_aux_basis, ri_hfx_basis, & 65 ri_hxc_basis, ri_xas_basis 66 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set 67 TYPE(qs_kind_type), POINTER :: qs_kind 68 69 IF (ncalls > 0) RETURN 70 ncalls = ncalls + 1 71 72 logger => cp_get_default_logger() 73 ounit = cp_logger_get_default_io_unit(logger) 74 75 CALL section_vals_val_get(base_section, "FILENAME", c_val=filename) 76 77 IF (ounit > 0) THEN 78 WRITE (UNIT=ounit, FMT='(/,(T2,A))') REPEAT("-", 79) 79 WRITE (UNIT=ounit, FMT='((T2,A,A))') "Print Basis Set File: ", TRIM(filename) 80 WRITE (UNIT=ounit, FMT='((T2,A))') REPEAT("-", 79) 81 CALL open_file(filename, unit_number=iunit, file_status="UNKNOWN", file_action="WRITE") 82 WRITE (UNIT=iunit, FMT="(A8,T11,A)") & 83 "# TITLE ", "Basis set file created by "//TRIM(cp2k_version)//" (revision "//TRIM(compile_revision)//")", & 84 "# AUTHOR", TRIM(r_user_name)//"@"//TRIM(r_host_name)//" "//r_datx(1:19) 85 86 ENDIF 87 88 CALL get_qs_env(qs_env, qs_kind_set=qs_kind_set, nkind=nkind) 89 DO ikind = 1, nkind 90 qs_kind => qs_kind_set(ikind) 91 CALL get_qs_kind(qs_kind, element_symbol=element_symbol) 92 NULLIFY (orb_basis, ri_aux_basis, lri_aux_basis, aux_fit_basis) 93 CALL get_qs_kind(qs_kind, basis_set=orb_basis, basis_type="ORB") 94 CALL get_qs_kind(qs_kind, basis_set=ri_aux_basis, basis_type="RI_AUX") 95 CALL get_qs_kind(qs_kind, basis_set=ri_hxc_basis, basis_type="RI_HXC") 96 CALL get_qs_kind(qs_kind, basis_set=ri_hfx_basis, basis_type="RI_HFX") 97 CALL get_qs_kind(qs_kind, basis_set=lri_aux_basis, basis_type="LRI_AUX") 98 CALL get_qs_kind(qs_kind, basis_set=aux_fit_basis, basis_type="AUX_FIT") 99 CALL get_qs_kind(qs_kind, basis_set=ri_xas_basis, basis_type="RI_XAS") 100 IF (ounit > 0) THEN 101 IF (ASSOCIATED(orb_basis)) THEN 102 bname = "local_orbital" 103 CALL basis_out(orb_basis, element_symbol, bname, iunit) 104 END IF 105 IF (ASSOCIATED(ri_aux_basis)) THEN 106 bname = "local_ri_aux" 107 CALL basis_out(ri_aux_basis, element_symbol, bname, iunit) 108 END IF 109 IF (ASSOCIATED(ri_hxc_basis)) THEN 110 bname = "local_ri_hxc" 111 CALL basis_out(ri_hxc_basis, element_symbol, bname, iunit) 112 END IF 113 IF (ASSOCIATED(lri_aux_basis)) THEN 114 bname = "local_lri_aux" 115 CALL basis_out(lri_aux_basis, element_symbol, bname, iunit) 116 END IF 117 IF (ASSOCIATED(aux_fit_basis)) THEN 118 bname = "local_aux_fit" 119 CALL basis_out(aux_fit_basis, element_symbol, bname, iunit) 120 END IF 121 IF (ASSOCIATED(ri_xas_basis)) THEN 122 bname = "local_ri_xas" 123 CALL basis_out(ri_xas_basis, element_symbol, bname, iunit) 124 END IF 125 IF (ASSOCIATED(ri_hfx_basis)) THEN 126 bname = "local_ri_hfx" 127 CALL basis_out(ri_hfx_basis, element_symbol, bname, iunit) 128 ENDIF 129 ENDIF 130 END DO 131 132 IF (ounit > 0) THEN 133 CALL close_file(iunit) 134 ENDIF 135 136 END SUBROUTINE print_basis_set_file 137 138! ************************************************************************************************** 139 140! ************************************************************************************************** 141!> \brief ... 142!> \param basis ... 143!> \param element_symbol ... 144!> \param bname ... 145!> \param iunit ... 146! ************************************************************************************************** 147 SUBROUTINE basis_out(basis, element_symbol, bname, iunit) 148 TYPE(gto_basis_set_type), POINTER :: basis 149 CHARACTER(LEN=*), INTENT(IN) :: element_symbol, bname 150 INTEGER, INTENT(IN) :: iunit 151 152 CHARACTER(len=*), PARAMETER :: routineN = 'basis_out', routineP = moduleN//':'//routineN 153 154 INTEGER :: ipgf, iset, ishell, ll, nset 155 INTEGER, DIMENSION(0:9) :: lset 156 INTEGER, DIMENSION(:), POINTER :: lmax, lmin, npgf, nshell 157 INTEGER, DIMENSION(:, :), POINTER :: l, n 158 REAL(KIND=dp), DIMENSION(:, :), POINTER :: zet 159 REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: gcc 160 161 WRITE (iunit, "(A1)") "#" 162 WRITE (iunit, "(A2,T5,A)") element_symbol, ADJUSTL(TRIM(bname)) 163 164 CALL get_gto_basis_set(basis, nset=nset, npgf=npgf, lmax=lmax, lmin=lmin, & 165 nshell=nshell, n=n, l=l, & 166 gcc=gcc, zet=zet) 167 168 WRITE (iunit, "(I5)") nset 169 DO iset = 1, nset 170 lset = 0 171 DO ishell = 1, nshell(iset) 172 ll = l(ishell, iset) 173 lset(ll) = lset(ll) + 1 174 END DO 175 WRITE (iunit, "(I5,2I3,I5,2X,10(I3))") n(1, iset), lmin(iset), lmax(iset), npgf(iset), & 176 (lset(ll), ll=lmin(iset), lmax(iset)) 177 DO ipgf = 1, npgf(iset) 178 WRITE (iunit, "(F20.10,50(F15.10))") zet(ipgf, iset), (gcc(ipgf, ishell, iset), ishell=1, nshell(iset)) 179 END DO 180 END DO 181 182 END SUBROUTINE basis_out 183 184! ************************************************************************************************** 185 186END MODULE basis_set_output 187