(git:374b731)
Loading...
Searching...
No Matches
qs_gcp_utils.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Set disperson types for DFT calculations
10!> \author JGH (04.2014)
11! **************************************************************************************************
13
25 USE kinds, ONLY: default_string_length,&
26 dp
27 USE mathconstants, ONLY: pi
32 USE qs_gcp_types, ONLY: qs_gcp_type
33 USE qs_kind_types, ONLY: get_qs_kind,&
35 USE sto_ng, ONLY: get_sto_ng
36#include "./base/base_uses.f90"
37
38 IMPLICIT NONE
39
40 PRIVATE
41
42 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_gcp_utils'
43
44 INTEGER, DIMENSION(106) :: nshell = (/ &
45 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, & ! 1-30
46 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 6, 6, 6, 6, 6, 6, & ! 31-60
47 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 7, & ! 61-90
48 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7/) ! 91-106
49
50 INTEGER, DIMENSION(106) :: nll = (/ &
51 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 2, & ! 1-30
52 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 4, 4, 4, 4, & ! 31-60
53 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 3, 3, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 4, 4, & ! 61-90
54 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 3, 3, 0, 0/) ! 91-106
55
56 !
57 ! Slater exponents for valence states
58 ! Aleksander Herman
59 ! Empirically adjusted and consistent set of EHT valence orbital parameters for all elements of the periodic table
60 ! Modelling and Simulation in Materials Science and Engineering, 12, 21-32 (2004)
61 !
62 ! Hydrogen uses 1.2000, not the original 1.0000
63 !
64 REAL(KIND=dp), DIMENSION(4, 106) :: sexp = reshape((/ &
65 1.2000, 0.0000, 0.0000, 0.0000, 1.6469, 0.0000, 0.0000, 0.0000, 0.6534, 0.5305, 0.0000, 0.0000, & ! 1 2 3
66 1.0365, 0.8994, 0.0000, 0.0000, 1.3990, 1.2685, 0.0000, 0.0000, 1.7210, 1.6105, 0.0000, 0.0000, & ! 4 5 6
67 2.0348, 1.9398, 0.0000, 0.0000, 2.2399, 2.0477, 0.0000, 0.0000, 2.5644, 2.4022, 0.0000, 0.0000, & ! 7 8 9
68 2.8812, 2.7421, 0.0000, 0.0000, 0.8675, 0.6148, 0.0000, 0.0000, 1.1935, 0.8809, 0.0000, 0.0000, & ! 10 11 12
69 1.5143, 1.1660, 0.0000, 0.0000, 1.7580, 1.4337, 0.0000, 0.0000, 1.9860, 1.6755, 0.0000, 0.0000, & ! 13 14 15
70 2.1362, 1.7721, 0.0000, 0.0000, 2.3617, 2.0176, 0.0000, 0.0000, 2.5796, 2.2501, 0.0000, 0.0000, & ! 16 17 18
71 0.9362, 0.6914, 0.0000, 0.0000, 1.2112, 0.9329, 0.0000, 0.0000, 1.2870, 0.9828, 2.4341, 0.0000, & ! 19 20 21
72 1.3416, 1.0104, 2.6439, 0.0000, 1.3570, 0.9947, 2.7809, 0.0000, 1.3804, 0.9784, 2.9775, 0.0000, & ! 22 23 24
73 1.4761, 1.0641, 3.2208, 0.0000, 1.5465, 1.1114, 3.4537, 0.0000, 1.5650, 1.1001, 3.6023, 0.0000, & ! 25 26 27
74 1.5532, 1.0594, 3.7017, 0.0000, 1.5791, 1.0527, 3.8962, 0.0000, 1.7778, 1.2448, 0.0000, 0.0000, & ! 28 29 30
75 2.0675, 1.5073, 0.0000, 0.0000, 2.2702, 1.7680, 0.0000, 0.0000, 2.4546, 1.9819, 0.0000, 0.0000, & ! 31 32 33
76 2.5680, 2.0548, 0.0000, 0.0000, 2.7523, 2.2652, 0.0000, 0.0000, 2.9299, 2.4617, 0.0000, 0.0000, & ! 34 35 36
77 1.0963, 0.7990, 0.0000, 0.0000, 1.3664, 1.0415, 0.0000, 0.0000, 1.4613, 1.1100, 2.1576, 0.0000, & ! 37 38 39
78 1.5393, 1.1647, 2.3831, 0.0000, 1.5926, 1.1738, 2.6256, 0.0000, 1.6579, 1.2186, 2.8241, 0.0000, & ! 40 41 42
79 1.6930, 1.2490, 2.9340, 0.0000, 1.7347, 1.2514, 3.1524, 0.0000, 1.7671, 1.2623, 3.3113, 0.0000, & ! 43 44 45
80 1.6261, 1.1221, 3.0858, 0.0000, 1.8184, 1.2719, 3.6171, 0.0000, 1.9900, 1.4596, 0.0000, 0.0000, & ! 46 47 48
81 2.4649, 1.6848, 0.0000, 0.0000, 2.4041, 1.9128, 0.0000, 0.0000, 2.5492, 2.0781, 0.0000, 0.0000, & ! 49 50 51
82 2.6576, 2.1718, 0.0000, 0.0000, 2.8080, 2.3390, 0.0000, 0.0000, 2.9595, 2.5074, 0.0000, 0.0000, & ! 52 53 54
83 1.1993, 0.8918, 0.0000, 0.0000, 1.4519, 1.1397, 0.0000, 0.0000, 1.5331, 1.1979, 2.2743, 4.4161, & ! 55 56 57
84 1.5379, 1.1930, 2.2912, 4.9478, 1.5162, 1.1834, 2.0558, 4.8982, 1.5322, 1.1923, 2.0718, 5.0744, & ! 58 59 60
85 1.5486, 1.2018, 2.0863, 5.2466, 1.5653, 1.2118, 2.0999, 5.4145, 1.5762, 1.2152, 2.0980, 5.5679, & ! 61 62 63
86 1.6703, 1.2874, 2.4862, 5.9888, 1.6186, 1.2460, 2.1383, 5.9040, 1.6358, 1.2570, 2.1472, 6.0598, & ! 64 65 66
87 1.6536, 1.2687, 2.1566, 6.2155, 1.6723, 1.2813, 2.1668, 6.3703, 1.6898, 1.2928, 2.1731, 6.5208, & ! 67 68 69
88 1.7063, 1.3030, 2.1754, 6.6686, 1.6647, 1.2167, 2.3795, 0.0000, 1.8411, 1.3822, 2.7702, 0.0000, & ! 70 71 72
89 1.9554, 1.4857, 3.0193, 0.0000, 2.0190, 1.5296, 3.1936, 0.0000, 2.0447, 1.5276, 3.3237, 0.0000, & ! 73 74 75
90 2.1361, 1.6102, 3.5241, 0.0000, 2.2167, 1.6814, 3.7077, 0.0000, 2.2646, 1.6759, 3.8996, 0.0000, & ! 76 77 78
91 2.3185, 1.7126, 4.0525, 0.0000, 2.4306, 1.8672, 0.0000, 0.0000, 2.5779, 1.9899, 0.0000, 0.0000, & ! 79 80 81
92 2.7241, 2.1837, 0.0000, 0.0000, 2.7869, 2.2146, 0.0000, 0.0000, 2.9312, 2.3830, 0.0000, 0.0000, & ! 82 83 84
93 3.1160, 2.6200, 0.0000, 0.0000, 3.2053, 2.6866, 0.0000, 0.0000, 1.4160, 1.0598, 0.0000, 0.0000, & ! 85 86 87
94 1.6336, 1.3011, 0.0000, 0.0000, 1.6540, 1.2890, 2.3740, 3.7960, 1.8381, 1.4726, 2.6584, 4.3613, & ! 88 89 90
95 1.7770, 1.4120, 2.5710, 4.5540, 1.8246, 1.4588, 2.6496, 4.7702, 1.8451, 1.4739, 2.6940, 4.9412, & ! 91 92 93
96 1.7983, 1.4366, 2.5123, 4.9882, 1.8011, 1.4317, 2.5170, 5.1301, 1.8408, 1.4418, 2.7349, 5.3476, & ! 94 95 96
97 1.8464, 1.4697, 2.5922, 5.4596, 1.8647, 1.4838, 2.6205, 5.6140, 1.8890, 1.5050, 2.6590, 5.7740, & ! 97 98 99
98 1.9070, 1.5190, 2.6850, 5.9220, 1.9240, 1.5320, 2.7090, 6.0690, 1.9400, 1.5440, 2.7300, 6.2130, & ! 100 101 102
99 2.1300, 1.7200, 2.9900, 0.0000, 1.9200, 1.4500, 2.9700, 0.0000, 0.0000, 0.0000, 0.0000, 0.0000, & ! 103 104 105
100 0.0000, 0.0000, 0.0000, 0.0000/), & ! 106
101 (/4, 106/))
102
103 PUBLIC :: qs_gcp_env_set, qs_gcp_init
104
105! **************************************************************************************************
106CONTAINS
107! **************************************************************************************************
108!> \brief ...
109!> \param gcp_env ...
110!> \param xc_section ...
111! **************************************************************************************************
112 SUBROUTINE qs_gcp_env_set(gcp_env, xc_section)
113 TYPE(qs_gcp_type), POINTER :: gcp_env
114 TYPE(section_vals_type), POINTER :: xc_section
115
116 CHARACTER(LEN=default_string_length), &
117 DIMENSION(:), POINTER :: tmpstringlist
118 INTEGER :: i_rep, n_rep
119 LOGICAL :: explicit
120 REAL(dp), POINTER :: params(:)
121 TYPE(section_vals_type), POINTER :: gcp_section
122
123 cpassert(ASSOCIATED(gcp_env))
124
125 gcp_section => section_vals_get_subs_vals(xc_section, "GCP_POTENTIAL")
126 CALL section_vals_get(gcp_section, explicit=explicit)
127 IF (explicit) THEN
128 CALL section_vals_val_get(gcp_section, "VERBOSE", l_val=gcp_env%verbose)
129 gcp_env%do_gcp = .true.
130 CALL section_vals_val_get(gcp_section, "PARAMETER_FILE_NAME", &
131 c_val=gcp_env%parameter_file_name)
132 CALL section_vals_val_get(gcp_section, "GLOBAL_PARAMETERS", r_vals=params)
133 gcp_env%sigma = params(1)
134 gcp_env%alpha = params(2)
135 gcp_env%beta = params(3)
136 gcp_env%eta = params(4)
137 ! eamiss definitions
138 CALL section_vals_val_get(gcp_section, "DELTA_ENERGY", n_rep_val=n_rep)
139 IF (n_rep > 0) THEN
140 ALLOCATE (gcp_env%kind_type(n_rep))
141 ALLOCATE (gcp_env%ea(n_rep))
142 DO i_rep = 1, n_rep
143 CALL section_vals_val_get(gcp_section, "DELTA_ENERGY", i_rep_val=i_rep, &
144 c_vals=tmpstringlist)
145 READ (tmpstringlist(1), *) gcp_env%kind_type(i_rep)
146 READ (tmpstringlist(2), *) gcp_env%ea(i_rep)
147 END DO
148 END IF
149 ELSE
150 gcp_env%do_gcp = .false.
151 END IF
152
153 END SUBROUTINE qs_gcp_env_set
154
155! **************************************************************************************************
156!> \brief ...
157!> \param qs_env ...
158!> \param gcp_env ...
159! **************************************************************************************************
160 SUBROUTINE qs_gcp_init(qs_env, gcp_env)
161 TYPE(qs_environment_type), POINTER :: qs_env
162 TYPE(qs_gcp_type), POINTER :: gcp_env
163
164 REAL(kind=dp), PARAMETER :: epsc = 1.e-6_dp
165
166 CHARACTER(LEN=10) :: aname
167 CHARACTER(LEN=2) :: element_symbol
168 INTEGER :: i, ikind, nbas, nel, nkind, nsto, za
169 LOGICAL :: at_end
170 REAL(kind=dp) :: ea
171 REAL(kind=dp), DIMENSION(10) :: al, cl
172 TYPE(gto_basis_set_type), POINTER :: orb_basis
173 TYPE(mp_para_env_type), POINTER :: para_env
174 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
175 TYPE(qs_kind_type), POINTER :: qs_kind
176
177 IF (gcp_env%do_gcp) THEN
178 CALL get_qs_env(qs_env, qs_kind_set=qs_kind_set, nkind=nkind)
179 ALLOCATE (gcp_env%gcp_kind(nkind))
180 DO ikind = 1, nkind
181 qs_kind => qs_kind_set(ikind)
182 gcp_env%gcp_kind(ikind)%rcsto = 0.0_dp
183 CALL get_qs_kind(qs_kind, element_symbol=element_symbol)
184 CALL get_ptable_info(element_symbol, number=za)
185 gcp_env%gcp_kind(ikind)%za = za
186 gcp_env%gcp_kind(ikind)%asto = gcp_env%eta*sum(sexp(1:4, za))/real(nll(za), kind=dp)
187 gcp_env%gcp_kind(ikind)%nq = nshell(za)
188 gcp_env%gcp_kind(ikind)%rcsto = ((nshell(za) - 1)*2.5_dp - log(epsc))/gcp_env%gcp_kind(ikind)%asto
189 ! basis
190 NULLIFY (orb_basis)
191 CALL get_qs_kind(qs_kind, basis_set=orb_basis, basis_type="ORB")
192 CALL get_gto_basis_set(gto_basis_set=orb_basis, nsgf=nbas)
193 nel = sum(qs_kind%elec_conf)
194 gcp_env%gcp_kind(ikind)%nbvirt = real(nbas, kind=dp) - 0.5_dp*real(nel, kind=dp)
195 ! STO-nG
196 nsto = SIZE(gcp_env%gcp_kind(ikind)%al)
197 CALL get_sto_ng(gcp_env%gcp_kind(ikind)%asto, nsto, nshell(za), 0, al, cl)
198 DO i = 1, nsto
199 gcp_env%gcp_kind(ikind)%al(i) = al(i)
200 gcp_env%gcp_kind(ikind)%cl(i) = cl(i)*(2._dp*al(i)/pi)**0.75_dp
201 END DO
202 END DO
203 ! eamiss from data file
204 IF (gcp_env%parameter_file_name /= "---") THEN
205 block
206 TYPE(cp_parser_type) :: parser
207 CALL get_qs_env(qs_env, para_env=para_env)
208 DO ikind = 1, nkind
209 qs_kind => qs_kind_set(ikind)
210 CALL get_qs_kind(qs_kind, element_symbol=element_symbol)
211 CALL get_ptable_info(element_symbol, number=za)
212 !
213 CALL parser_create(parser, gcp_env%parameter_file_name, para_env=para_env)
214 ea = 0.0_dp
215 DO
216 at_end = .false.
217 CALL parser_get_next_line(parser, 1, at_end)
218 IF (at_end) EXIT
219 CALL parser_get_object(parser, aname)
220 IF (trim(aname) == element_symbol) THEN
221 CALL parser_get_object(parser, ea)
222 EXIT
223 END IF
224 END DO
225 CALL parser_release(parser)
226 gcp_env%gcp_kind(ikind)%eamiss = ea
227 END DO
228 END block
229 END IF
230 !
231 ! eamiss from input
232 IF (ASSOCIATED(gcp_env%kind_type)) THEN
233 DO i = 1, SIZE(gcp_env%kind_type)
234 IF (trim(gcp_env%kind_type(i)) == "XX") cycle
235 element_symbol = trim(gcp_env%kind_type(i))
236 CALL get_ptable_info(element_symbol, number=za)
237 ea = gcp_env%ea(i)
238 DO ikind = 1, nkind
239 IF (za == gcp_env%gcp_kind(ikind)%za) THEN
240 gcp_env%gcp_kind(ikind)%eamiss = ea
241 END IF
242 END DO
243 END DO
244 END IF
245 END IF
246
247 END SUBROUTINE qs_gcp_init
248! **************************************************************************************************
249
250END MODULE qs_gcp_utils
251
subroutine, public get_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, last_cgf, last_sgf, n, gcc, maxco, maxl, maxpgf, maxsgf_set, maxshell, maxso, nco_sum, npgf_sum, nshell_sum, maxder, short_kind_radius)
...
Utility routines to read data from files. Kept as close as possible to the old parser because.
subroutine, public parser_get_next_line(parser, nline, at_end)
Read the next input line and broadcast the input information. Skip (nline-1) lines and skip also all ...
Utility routines to read data from files. Kept as close as possible to the old parser because.
subroutine, public parser_release(parser)
releases the parser
subroutine, public parser_create(parser, file_name, unit_nr, para_env, end_section_label, separator_chars, comment_char, continuation_char, quote_char, section_char, parse_white_lines, initial_variables, apply_preprocessing)
Start a parser run. Initial variables allow to @SET stuff before opening the file.
objects that represent the structure of input sections and the data contained in an input section
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
integer, parameter, public default_string_length
Definition kinds.F:57
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
Interface to the message passing library MPI.
Periodic Table related data definitions.
subroutine, public get_ptable_info(symbol, number, amass, ielement, covalent_radius, metallic_radius, vdw_radius, found)
Pass information about the kind given the element symbol.
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_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_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, ecoul_1c, rho0_s_rs, rho0_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Get the QUICKSTEP environment.
Definition of gCP types for DFT calculations.
Set disperson types for DFT calculations.
subroutine, public qs_gcp_env_set(gcp_env, xc_section)
...
subroutine, public qs_gcp_init(qs_env, gcp_env)
...
Define the quickstep kind type and their sub types.
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_r3d_rs_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, u_of_dft_plus_u, j_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, j0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, floating, name, element_symbol, pao_basis_size, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
subroutine, public get_sto_ng(zeta, n, nq, lq, alpha, coef)
return STO-NG parameters; INPUT: zeta (Slater exponent) n (Expansion length) nq (principle quantum nu...
Definition sto_ng.F:54
stores all the informations relevant to an mpi environment
Provides all information about a quickstep kind.