Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 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 : ! **************************************************************************************************
12 : MODULE qs_gcp_utils
13 :
14 : USE basis_set_types, ONLY: get_gto_basis_set,&
15 : gto_basis_set_type
16 : USE cp_parser_methods, ONLY: parser_get_next_line,&
17 : parser_get_object
18 : USE cp_parser_types, ONLY: cp_parser_type,&
19 : parser_create,&
20 : parser_release
21 : USE input_section_types, ONLY: section_vals_get,&
22 : section_vals_get_subs_vals,&
23 : section_vals_type,&
24 : section_vals_val_get
25 : USE kinds, ONLY: default_string_length,&
26 : dp
27 : USE mathconstants, ONLY: pi
28 : USE message_passing, ONLY: mp_para_env_type
29 : USE periodic_table, ONLY: get_ptable_info
30 : USE qs_environment_types, ONLY: get_qs_env,&
31 : qs_environment_type
32 : USE qs_gcp_types, ONLY: qs_gcp_type
33 : USE qs_kind_types, ONLY: get_qs_kind,&
34 : qs_kind_type
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), PARAMETER, PRIVATE :: sexp = RESHAPE([ &
65 : 1.2000_dp, 0.0000_dp, 0.0000_dp, 0.0000_dp, 1.6469_dp, 0.0000_dp, 0.0000_dp, 0.0000_dp, & ! 1 2
66 : 0.6534_dp, 0.5305_dp, 0.0000_dp, 0.0000_dp, 1.0365_dp, 0.8994_dp, 0.0000_dp, 0.0000_dp, & ! 3 4
67 : 1.3990_dp, 1.2685_dp, 0.0000_dp, 0.0000_dp, 1.7210_dp, 1.6105_dp, 0.0000_dp, 0.0000_dp, & ! 5 6
68 : 2.0348_dp, 1.9398_dp, 0.0000_dp, 0.0000_dp, 2.2399_dp, 2.0477_dp, 0.0000_dp, 0.0000_dp, & ! 7 8
69 : 2.5644_dp, 2.4022_dp, 0.0000_dp, 0.0000_dp, 2.8812_dp, 2.7421_dp, 0.0000_dp, 0.0000_dp, & ! 9 10
70 : 0.8675_dp, 0.6148_dp, 0.0000_dp, 0.0000_dp, 1.1935_dp, 0.8809_dp, 0.0000_dp, 0.0000_dp, & ! 11 12
71 : 1.5143_dp, 1.1660_dp, 0.0000_dp, 0.0000_dp, 1.7580_dp, 1.4337_dp, 0.0000_dp, 0.0000_dp, & ! 13 14
72 : 1.9860_dp, 1.6755_dp, 0.0000_dp, 0.0000_dp, 2.1362_dp, 1.7721_dp, 0.0000_dp, 0.0000_dp, & ! 15 16
73 : 2.3617_dp, 2.0176_dp, 0.0000_dp, 0.0000_dp, 2.5796_dp, 2.2501_dp, 0.0000_dp, 0.0000_dp, & ! 17 18
74 : 0.9362_dp, 0.6914_dp, 0.0000_dp, 0.0000_dp, 1.2112_dp, 0.9329_dp, 0.0000_dp, 0.0000_dp, & ! 19 20
75 : 1.2870_dp, 0.9828_dp, 2.4341_dp, 0.0000_dp, 1.3416_dp, 1.0104_dp, 2.6439_dp, 0.0000_dp, & ! 21 22
76 : 1.3570_dp, 0.9947_dp, 2.7809_dp, 0.0000_dp, 1.3804_dp, 0.9784_dp, 2.9775_dp, 0.0000_dp, & ! 23 24
77 : 1.4761_dp, 1.0641_dp, 3.2208_dp, 0.0000_dp, 1.5465_dp, 1.1114_dp, 3.4537_dp, 0.0000_dp, & ! 25 26
78 : 1.5650_dp, 1.1001_dp, 3.6023_dp, 0.0000_dp, 1.5532_dp, 1.0594_dp, 3.7017_dp, 0.0000_dp, & ! 27 28
79 : 1.5791_dp, 1.0527_dp, 3.8962_dp, 0.0000_dp, 1.7778_dp, 1.2448_dp, 0.0000_dp, 0.0000_dp, & ! 29 30
80 : 2.0675_dp, 1.5073_dp, 0.0000_dp, 0.0000_dp, 2.2702_dp, 1.7680_dp, 0.0000_dp, 0.0000_dp, & ! 31 32
81 : 2.4546_dp, 1.9819_dp, 0.0000_dp, 0.0000_dp, 2.5680_dp, 2.0548_dp, 0.0000_dp, 0.0000_dp, & ! 33 34
82 : 2.7523_dp, 2.2652_dp, 0.0000_dp, 0.0000_dp, 2.9299_dp, 2.4617_dp, 0.0000_dp, 0.0000_dp, & ! 35 36
83 : 1.0963_dp, 0.7990_dp, 0.0000_dp, 0.0000_dp, 1.3664_dp, 1.0415_dp, 0.0000_dp, 0.0000_dp, & ! 37 38
84 : 1.4613_dp, 1.1100_dp, 2.1576_dp, 0.0000_dp, 1.5393_dp, 1.1647_dp, 2.3831_dp, 0.0000_dp, & ! 39 40
85 : 1.5926_dp, 1.1738_dp, 2.6256_dp, 0.0000_dp, 1.6579_dp, 1.2186_dp, 2.8241_dp, 0.0000_dp, & ! 41 42
86 : 1.6930_dp, 1.2490_dp, 2.9340_dp, 0.0000_dp, 1.7347_dp, 1.2514_dp, 3.1524_dp, 0.0000_dp, & ! 43 44
87 : 1.7671_dp, 1.2623_dp, 3.3113_dp, 0.0000_dp, 1.6261_dp, 1.1221_dp, 3.0858_dp, 0.0000_dp, & ! 45 46
88 : 1.8184_dp, 1.2719_dp, 3.6171_dp, 0.0000_dp, 1.9900_dp, 1.4596_dp, 0.0000_dp, 0.0000_dp, & ! 47 48
89 : 2.4649_dp, 1.6848_dp, 0.0000_dp, 0.0000_dp, 2.4041_dp, 1.9128_dp, 0.0000_dp, 0.0000_dp, & ! 49 50
90 : 2.5492_dp, 2.0781_dp, 0.0000_dp, 0.0000_dp, 2.6576_dp, 2.1718_dp, 0.0000_dp, 0.0000_dp, & ! 51 52
91 : 2.8080_dp, 2.3390_dp, 0.0000_dp, 0.0000_dp, 2.9595_dp, 2.5074_dp, 0.0000_dp, 0.0000_dp, & ! 53 54
92 : 1.1993_dp, 0.8918_dp, 0.0000_dp, 0.0000_dp, 1.4519_dp, 1.1397_dp, 0.0000_dp, 0.0000_dp, & ! 55 56
93 : 1.5331_dp, 1.1979_dp, 2.2743_dp, 4.4161_dp, 1.5379_dp, 1.1930_dp, 2.2912_dp, 4.9478_dp, & ! 57 58
94 : 1.5162_dp, 1.1834_dp, 2.0558_dp, 4.8982_dp, 1.5322_dp, 1.1923_dp, 2.0718_dp, 5.0744_dp, & ! 59 60
95 : 1.5486_dp, 1.2018_dp, 2.0863_dp, 5.2466_dp, 1.5653_dp, 1.2118_dp, 2.0999_dp, 5.4145_dp, & ! 61 62
96 : 1.5762_dp, 1.2152_dp, 2.0980_dp, 5.5679_dp, 1.6703_dp, 1.2874_dp, 2.4862_dp, 5.9888_dp, & ! 63 64
97 : 1.6186_dp, 1.2460_dp, 2.1383_dp, 5.9040_dp, 1.6358_dp, 1.2570_dp, 2.1472_dp, 6.0598_dp, & ! 65 66
98 : 1.6536_dp, 1.2687_dp, 2.1566_dp, 6.2155_dp, 1.6723_dp, 1.2813_dp, 2.1668_dp, 6.3703_dp, & ! 67 68
99 : 1.6898_dp, 1.2928_dp, 2.1731_dp, 6.5208_dp, 1.7063_dp, 1.3030_dp, 2.1754_dp, 6.6686_dp, & ! 69 70
100 : 1.6647_dp, 1.2167_dp, 2.3795_dp, 0.0000_dp, 1.8411_dp, 1.3822_dp, 2.7702_dp, 0.0000_dp, & ! 71 72
101 : 1.9554_dp, 1.4857_dp, 3.0193_dp, 0.0000_dp, 2.0190_dp, 1.5296_dp, 3.1936_dp, 0.0000_dp, & ! 73 74
102 : 2.0447_dp, 1.5276_dp, 3.3237_dp, 0.0000_dp, 2.1361_dp, 1.6102_dp, 3.5241_dp, 0.0000_dp, & ! 75 76
103 : 2.2167_dp, 1.6814_dp, 3.7077_dp, 0.0000_dp, 2.2646_dp, 1.6759_dp, 3.8996_dp, 0.0000_dp, & ! 77 78
104 : 2.3185_dp, 1.7126_dp, 4.0525_dp, 0.0000_dp, 2.4306_dp, 1.8672_dp, 0.0000_dp, 0.0000_dp, & ! 79 80
105 : 2.5779_dp, 1.9899_dp, 0.0000_dp, 0.0000_dp, 2.7241_dp, 2.1837_dp, 0.0000_dp, 0.0000_dp, & ! 81 82
106 : 2.7869_dp, 2.2146_dp, 0.0000_dp, 0.0000_dp, 2.9312_dp, 2.3830_dp, 0.0000_dp, 0.0000_dp, & ! 83 84
107 : 3.1160_dp, 2.6200_dp, 0.0000_dp, 0.0000_dp, 3.2053_dp, 2.6866_dp, 0.0000_dp, 0.0000_dp, & ! 85 86
108 : 1.4160_dp, 1.0598_dp, 0.0000_dp, 0.0000_dp, 1.6336_dp, 1.3011_dp, 0.0000_dp, 0.0000_dp, & ! 87 88
109 : 1.6540_dp, 1.2890_dp, 2.3740_dp, 3.7960_dp, 1.8381_dp, 1.4726_dp, 2.6584_dp, 4.3613_dp, & ! 89 90
110 : 1.7770_dp, 1.4120_dp, 2.5710_dp, 4.5540_dp, 1.8246_dp, 1.4588_dp, 2.6496_dp, 4.7702_dp, & ! 91 92
111 : 1.8451_dp, 1.4739_dp, 2.6940_dp, 4.9412_dp, 1.7983_dp, 1.4366_dp, 2.5123_dp, 4.9882_dp, & ! 93 94
112 : 1.8011_dp, 1.4317_dp, 2.5170_dp, 5.1301_dp, 1.8408_dp, 1.4418_dp, 2.7349_dp, 5.3476_dp, & ! 95 96
113 : 1.8464_dp, 1.4697_dp, 2.5922_dp, 5.4596_dp, 1.8647_dp, 1.4838_dp, 2.6205_dp, 5.6140_dp, & ! 97 98
114 : 1.8890_dp, 1.5050_dp, 2.6590_dp, 5.7740_dp, 1.9070_dp, 1.5190_dp, 2.6850_dp, 5.9220_dp, & ! 99 100
115 : 1.9240_dp, 1.5320_dp, 2.7090_dp, 6.0690_dp, 1.9400_dp, 1.5440_dp, 2.7300_dp, 6.2130_dp, & ! 101 102
116 : 2.1300_dp, 1.7200_dp, 2.9900_dp, 0.0000_dp, 1.9200_dp, 1.4500_dp, 2.9700_dp, 0.0000_dp, & ! 103 104
117 : 0.0000_dp, 0.0000_dp, 0.0000_dp, 0.0000_dp, 0.0000_dp, 0.0000_dp, 0.0000_dp, 0.0000_dp], & ! 105 106
118 : [4, 106])
119 :
120 : PUBLIC :: qs_gcp_env_set, qs_gcp_init
121 :
122 : ! **************************************************************************************************
123 : CONTAINS
124 : ! **************************************************************************************************
125 : !> \brief ...
126 : !> \param gcp_env ...
127 : !> \param xc_section ...
128 : ! **************************************************************************************************
129 13268 : SUBROUTINE qs_gcp_env_set(gcp_env, xc_section)
130 : TYPE(qs_gcp_type), POINTER :: gcp_env
131 : TYPE(section_vals_type), POINTER :: xc_section
132 :
133 : CHARACTER(LEN=default_string_length), &
134 6634 : DIMENSION(:), POINTER :: tmpstringlist
135 : INTEGER :: i_rep, n_rep
136 : LOGICAL :: explicit
137 6634 : REAL(dp), POINTER :: params(:)
138 : TYPE(section_vals_type), POINTER :: gcp_section
139 :
140 0 : CPASSERT(ASSOCIATED(gcp_env))
141 :
142 6634 : gcp_section => section_vals_get_subs_vals(xc_section, "GCP_POTENTIAL")
143 6634 : CALL section_vals_get(gcp_section, explicit=explicit)
144 6634 : IF (explicit) THEN
145 6 : CALL section_vals_val_get(gcp_section, "VERBOSE", l_val=gcp_env%verbose)
146 6 : gcp_env%do_gcp = .TRUE.
147 : CALL section_vals_val_get(gcp_section, "PARAMETER_FILE_NAME", &
148 6 : c_val=gcp_env%parameter_file_name)
149 6 : CALL section_vals_val_get(gcp_section, "GLOBAL_PARAMETERS", r_vals=params)
150 6 : gcp_env%sigma = params(1)
151 6 : gcp_env%alpha = params(2)
152 6 : gcp_env%beta = params(3)
153 6 : gcp_env%eta = params(4)
154 : ! eamiss definitions
155 6 : CALL section_vals_val_get(gcp_section, "DELTA_ENERGY", n_rep_val=n_rep)
156 6 : IF (n_rep > 0) THEN
157 18 : ALLOCATE (gcp_env%kind_type(n_rep))
158 18 : ALLOCATE (gcp_env%ea(n_rep))
159 26 : DO i_rep = 1, n_rep
160 : CALL section_vals_val_get(gcp_section, "DELTA_ENERGY", i_rep_val=i_rep, &
161 20 : c_vals=tmpstringlist)
162 20 : READ (tmpstringlist(1), *) gcp_env%kind_type(i_rep)
163 26 : READ (tmpstringlist(2), *) gcp_env%ea(i_rep)
164 : END DO
165 : END IF
166 : ELSE
167 6628 : gcp_env%do_gcp = .FALSE.
168 : END IF
169 :
170 6634 : END SUBROUTINE qs_gcp_env_set
171 :
172 : ! **************************************************************************************************
173 : !> \brief ...
174 : !> \param qs_env ...
175 : !> \param gcp_env ...
176 : ! **************************************************************************************************
177 6634 : SUBROUTINE qs_gcp_init(qs_env, gcp_env)
178 : TYPE(qs_environment_type), POINTER :: qs_env
179 : TYPE(qs_gcp_type), POINTER :: gcp_env
180 :
181 : REAL(KIND=dp), PARAMETER :: epsc = 1.e-6_dp
182 :
183 : CHARACTER(LEN=10) :: aname
184 : CHARACTER(LEN=2) :: element_symbol
185 : INTEGER :: i, ikind, nbas, nel, nkind, nsto, za
186 : LOGICAL :: at_end
187 : REAL(KIND=dp) :: ea
188 : REAL(KIND=dp), DIMENSION(10) :: al, cl
189 : TYPE(gto_basis_set_type), POINTER :: orb_basis
190 : TYPE(mp_para_env_type), POINTER :: para_env
191 6634 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
192 : TYPE(qs_kind_type), POINTER :: qs_kind
193 :
194 6634 : IF (gcp_env%do_gcp) THEN
195 6 : CALL get_qs_env(qs_env, qs_kind_set=qs_kind_set, nkind=nkind)
196 98 : ALLOCATE (gcp_env%gcp_kind(nkind))
197 14 : DO ikind = 1, nkind
198 8 : qs_kind => qs_kind_set(ikind)
199 8 : gcp_env%gcp_kind(ikind)%rcsto = 0.0_dp
200 8 : CALL get_qs_kind(qs_kind, element_symbol=element_symbol)
201 8 : CALL get_ptable_info(element_symbol, number=za)
202 8 : gcp_env%gcp_kind(ikind)%za = za
203 40 : gcp_env%gcp_kind(ikind)%asto = gcp_env%eta*SUM(sexp(1:4, za))/REAL(nll(za), KIND=dp)
204 8 : gcp_env%gcp_kind(ikind)%nq = nshell(za)
205 8 : gcp_env%gcp_kind(ikind)%rcsto = ((nshell(za) - 1)*2.5_dp - LOG(epsc))/gcp_env%gcp_kind(ikind)%asto
206 : ! basis
207 8 : NULLIFY (orb_basis)
208 8 : CALL get_qs_kind(qs_kind, basis_set=orb_basis, basis_type="ORB")
209 8 : CALL get_gto_basis_set(gto_basis_set=orb_basis, nsgf=nbas)
210 32 : nel = SUM(qs_kind%elec_conf)
211 8 : gcp_env%gcp_kind(ikind)%nbvirt = REAL(nbas, KIND=dp) - 0.5_dp*REAL(nel, KIND=dp)
212 : ! STO-nG
213 8 : nsto = SIZE(gcp_env%gcp_kind(ikind)%al)
214 8 : CALL get_sto_ng(gcp_env%gcp_kind(ikind)%asto, nsto, nshell(za), 0, al, cl)
215 78 : DO i = 1, nsto
216 48 : gcp_env%gcp_kind(ikind)%al(i) = al(i)
217 56 : gcp_env%gcp_kind(ikind)%cl(i) = cl(i)*(2._dp*al(i)/pi)**0.75_dp
218 : END DO
219 : END DO
220 : ! eamiss from data file
221 6 : IF (gcp_env%parameter_file_name /= "---") THEN
222 : BLOCK
223 : TYPE(cp_parser_type) :: parser
224 0 : CALL get_qs_env(qs_env, para_env=para_env)
225 0 : DO ikind = 1, nkind
226 0 : qs_kind => qs_kind_set(ikind)
227 0 : CALL get_qs_kind(qs_kind, element_symbol=element_symbol)
228 0 : CALL get_ptable_info(element_symbol, number=za)
229 : !
230 0 : CALL parser_create(parser, gcp_env%parameter_file_name, para_env=para_env)
231 0 : ea = 0.0_dp
232 : DO
233 : at_end = .FALSE.
234 0 : CALL parser_get_next_line(parser, 1, at_end)
235 0 : IF (at_end) EXIT
236 0 : CALL parser_get_object(parser, aname)
237 0 : IF (TRIM(aname) == element_symbol) THEN
238 0 : CALL parser_get_object(parser, ea)
239 0 : EXIT
240 : END IF
241 : END DO
242 0 : CALL parser_release(parser)
243 0 : gcp_env%gcp_kind(ikind)%eamiss = ea
244 : END DO
245 : END BLOCK
246 : END IF
247 : !
248 : ! eamiss from input
249 6 : IF (ASSOCIATED(gcp_env%kind_type)) THEN
250 26 : DO i = 1, SIZE(gcp_env%kind_type)
251 20 : IF (TRIM(gcp_env%kind_type(i)) == "XX") CYCLE
252 20 : element_symbol = TRIM(gcp_env%kind_type(i))
253 20 : CALL get_ptable_info(element_symbol, number=za)
254 20 : ea = gcp_env%ea(i)
255 50 : DO ikind = 1, nkind
256 44 : IF (za == gcp_env%gcp_kind(ikind)%za) THEN
257 8 : gcp_env%gcp_kind(ikind)%eamiss = ea
258 : END IF
259 : END DO
260 : END DO
261 : END IF
262 : END IF
263 :
264 6634 : END SUBROUTINE qs_gcp_init
265 : ! **************************************************************************************************
266 :
267 : END MODULE qs_gcp_utils
268 :
|