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 : !> \par History
10 : !> Add CP2K error reporting, new add_force routine [07.2014,JGH]
11 : !> \author MK (03.06.2002)
12 : ! **************************************************************************************************
13 : MODULE qs_force_types
14 :
15 : USE atomic_kind_types, ONLY: atomic_kind_type,&
16 : get_atomic_kind
17 : USE cp_log_handling, ONLY: cp_get_default_logger,&
18 : cp_logger_get_default_io_unit,&
19 : cp_logger_type
20 : USE kinds, ONLY: dp
21 : USE message_passing, ONLY: mp_para_env_type
22 : #include "./base/base_uses.f90"
23 :
24 : IMPLICIT NONE
25 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_force_types'
26 : PRIVATE
27 :
28 : TYPE qs_force_type
29 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: all_potential => NULL(), &
30 : cneo_potential => NULL(), &
31 : core_overlap => NULL(), &
32 : gth_ppl => NULL(), &
33 : gth_nlcc => NULL(), &
34 : gth_ppnl => NULL(), &
35 : kinetic => NULL(), &
36 : overlap => NULL(), &
37 : overlap_admm => NULL(), &
38 : rho_core => NULL(), &
39 : rho_elec => NULL(), &
40 : rho_lri_elec => NULL(), &
41 : rho_cneo_nuc => NULL(), &
42 : vhxc_atom => NULL(), &
43 : g0s_Vh_elec => NULL(), &
44 : repulsive => NULL(), &
45 : dispersion => NULL(), &
46 : gcp => NULL(), &
47 : other => NULL(), &
48 : ch_pulay => NULL(), &
49 : fock_4c => NULL(), &
50 : ehrenfest => NULL(), &
51 : efield => NULL(), &
52 : eev => NULL(), &
53 : mp2_non_sep => NULL(), &
54 : total => NULL()
55 : END TYPE qs_force_type
56 :
57 : PUBLIC :: qs_force_type
58 :
59 : PUBLIC :: allocate_qs_force, &
60 : add_qs_force, &
61 : deallocate_qs_force, &
62 : replicate_qs_force, &
63 : sum_qs_force, &
64 : get_qs_force, &
65 : put_qs_force, &
66 : total_qs_force, &
67 : zero_qs_force, &
68 : write_forces_debug
69 :
70 : CONTAINS
71 :
72 : ! **************************************************************************************************
73 : !> \brief Allocate a Quickstep force data structure.
74 : !> \param qs_force ...
75 : !> \param natom_of_kind ...
76 : !> \date 05.06.2002
77 : !> \author MK
78 : !> \version 1.0
79 : ! **************************************************************************************************
80 5051 : SUBROUTINE allocate_qs_force(qs_force, natom_of_kind)
81 :
82 : TYPE(qs_force_type), DIMENSION(:), POINTER :: qs_force
83 : INTEGER, DIMENSION(:), INTENT(IN) :: natom_of_kind
84 :
85 : INTEGER :: ikind, n, nkind
86 :
87 5051 : IF (ASSOCIATED(qs_force)) CALL deallocate_qs_force(qs_force)
88 :
89 5051 : nkind = SIZE(natom_of_kind)
90 :
91 24869 : ALLOCATE (qs_force(nkind))
92 :
93 14767 : DO ikind = 1, nkind
94 9716 : n = natom_of_kind(ikind)
95 29148 : ALLOCATE (qs_force(ikind)%all_potential(3, n))
96 19432 : ALLOCATE (qs_force(ikind)%cneo_potential(3, n))
97 19432 : ALLOCATE (qs_force(ikind)%core_overlap(3, n))
98 19432 : ALLOCATE (qs_force(ikind)%gth_ppl(3, n))
99 19432 : ALLOCATE (qs_force(ikind)%gth_nlcc(3, n))
100 19432 : ALLOCATE (qs_force(ikind)%gth_ppnl(3, n))
101 19432 : ALLOCATE (qs_force(ikind)%kinetic(3, n))
102 19432 : ALLOCATE (qs_force(ikind)%overlap(3, n))
103 19432 : ALLOCATE (qs_force(ikind)%overlap_admm(3, n))
104 19432 : ALLOCATE (qs_force(ikind)%rho_core(3, n))
105 19432 : ALLOCATE (qs_force(ikind)%rho_elec(3, n))
106 19432 : ALLOCATE (qs_force(ikind)%rho_lri_elec(3, n))
107 19432 : ALLOCATE (qs_force(ikind)%rho_cneo_nuc(3, n))
108 19432 : ALLOCATE (qs_force(ikind)%vhxc_atom(3, n))
109 19432 : ALLOCATE (qs_force(ikind)%g0s_Vh_elec(3, n))
110 19432 : ALLOCATE (qs_force(ikind)%repulsive(3, n))
111 19432 : ALLOCATE (qs_force(ikind)%dispersion(3, n))
112 19432 : ALLOCATE (qs_force(ikind)%gcp(3, n))
113 19432 : ALLOCATE (qs_force(ikind)%other(3, n))
114 19432 : ALLOCATE (qs_force(ikind)%ch_pulay(3, n))
115 19432 : ALLOCATE (qs_force(ikind)%ehrenfest(3, n))
116 19432 : ALLOCATE (qs_force(ikind)%efield(3, n))
117 19432 : ALLOCATE (qs_force(ikind)%eev(3, n))
118 : ! Always initialize ch_pulay to zero..
119 113056 : qs_force(ikind)%ch_pulay = 0.0_dp
120 19432 : ALLOCATE (qs_force(ikind)%fock_4c(3, n))
121 19432 : ALLOCATE (qs_force(ikind)%mp2_non_sep(3, n))
122 24483 : ALLOCATE (qs_force(ikind)%total(3, n))
123 : END DO
124 :
125 5051 : END SUBROUTINE allocate_qs_force
126 :
127 : ! **************************************************************************************************
128 : !> \brief Deallocate a Quickstep force data structure.
129 : !> \param qs_force ...
130 : !> \date 05.06.2002
131 : !> \author MK
132 : !> \version 1.0
133 : ! **************************************************************************************************
134 5051 : SUBROUTINE deallocate_qs_force(qs_force)
135 :
136 : TYPE(qs_force_type), DIMENSION(:), POINTER :: qs_force
137 :
138 : INTEGER :: ikind, nkind
139 :
140 5051 : CPASSERT(ASSOCIATED(qs_force))
141 :
142 5051 : nkind = SIZE(qs_force)
143 :
144 14767 : DO ikind = 1, nkind
145 :
146 9716 : IF (ASSOCIATED(qs_force(ikind)%all_potential)) THEN
147 9716 : DEALLOCATE (qs_force(ikind)%all_potential)
148 : END IF
149 :
150 9716 : IF (ASSOCIATED(qs_force(ikind)%cneo_potential)) THEN
151 9716 : DEALLOCATE (qs_force(ikind)%cneo_potential)
152 : END IF
153 :
154 9716 : IF (ASSOCIATED(qs_force(ikind)%core_overlap)) THEN
155 9716 : DEALLOCATE (qs_force(ikind)%core_overlap)
156 : END IF
157 :
158 9716 : IF (ASSOCIATED(qs_force(ikind)%gth_ppl)) THEN
159 9716 : DEALLOCATE (qs_force(ikind)%gth_ppl)
160 : END IF
161 :
162 9716 : IF (ASSOCIATED(qs_force(ikind)%gth_nlcc)) THEN
163 9716 : DEALLOCATE (qs_force(ikind)%gth_nlcc)
164 : END IF
165 :
166 9716 : IF (ASSOCIATED(qs_force(ikind)%gth_ppnl)) THEN
167 9716 : DEALLOCATE (qs_force(ikind)%gth_ppnl)
168 : END IF
169 :
170 9716 : IF (ASSOCIATED(qs_force(ikind)%kinetic)) THEN
171 9716 : DEALLOCATE (qs_force(ikind)%kinetic)
172 : END IF
173 :
174 9716 : IF (ASSOCIATED(qs_force(ikind)%overlap)) THEN
175 9716 : DEALLOCATE (qs_force(ikind)%overlap)
176 : END IF
177 :
178 9716 : IF (ASSOCIATED(qs_force(ikind)%overlap_admm)) THEN
179 9716 : DEALLOCATE (qs_force(ikind)%overlap_admm)
180 : END IF
181 :
182 9716 : IF (ASSOCIATED(qs_force(ikind)%rho_core)) THEN
183 9716 : DEALLOCATE (qs_force(ikind)%rho_core)
184 : END IF
185 :
186 9716 : IF (ASSOCIATED(qs_force(ikind)%rho_elec)) THEN
187 9716 : DEALLOCATE (qs_force(ikind)%rho_elec)
188 : END IF
189 9716 : IF (ASSOCIATED(qs_force(ikind)%rho_lri_elec)) THEN
190 9716 : DEALLOCATE (qs_force(ikind)%rho_lri_elec)
191 : END IF
192 :
193 9716 : IF (ASSOCIATED(qs_force(ikind)%rho_cneo_nuc)) THEN
194 9716 : DEALLOCATE (qs_force(ikind)%rho_cneo_nuc)
195 : END IF
196 :
197 9716 : IF (ASSOCIATED(qs_force(ikind)%vhxc_atom)) THEN
198 9716 : DEALLOCATE (qs_force(ikind)%vhxc_atom)
199 : END IF
200 :
201 9716 : IF (ASSOCIATED(qs_force(ikind)%g0s_Vh_elec)) THEN
202 9716 : DEALLOCATE (qs_force(ikind)%g0s_Vh_elec)
203 : END IF
204 :
205 9716 : IF (ASSOCIATED(qs_force(ikind)%repulsive)) THEN
206 9716 : DEALLOCATE (qs_force(ikind)%repulsive)
207 : END IF
208 :
209 9716 : IF (ASSOCIATED(qs_force(ikind)%dispersion)) THEN
210 9716 : DEALLOCATE (qs_force(ikind)%dispersion)
211 : END IF
212 :
213 9716 : IF (ASSOCIATED(qs_force(ikind)%gcp)) THEN
214 9716 : DEALLOCATE (qs_force(ikind)%gcp)
215 : END IF
216 :
217 9716 : IF (ASSOCIATED(qs_force(ikind)%other)) THEN
218 9716 : DEALLOCATE (qs_force(ikind)%other)
219 : END IF
220 :
221 9716 : IF (ASSOCIATED(qs_force(ikind)%total)) THEN
222 9716 : DEALLOCATE (qs_force(ikind)%total)
223 : END IF
224 :
225 9716 : IF (ASSOCIATED(qs_force(ikind)%ch_pulay)) THEN
226 9716 : DEALLOCATE (qs_force(ikind)%ch_pulay)
227 : END IF
228 :
229 9716 : IF (ASSOCIATED(qs_force(ikind)%fock_4c)) THEN
230 9716 : DEALLOCATE (qs_force(ikind)%fock_4c)
231 : END IF
232 :
233 9716 : IF (ASSOCIATED(qs_force(ikind)%mp2_non_sep)) THEN
234 9716 : DEALLOCATE (qs_force(ikind)%mp2_non_sep)
235 : END IF
236 :
237 9716 : IF (ASSOCIATED(qs_force(ikind)%ehrenfest)) THEN
238 9716 : DEALLOCATE (qs_force(ikind)%ehrenfest)
239 : END IF
240 :
241 9716 : IF (ASSOCIATED(qs_force(ikind)%efield)) THEN
242 9716 : DEALLOCATE (qs_force(ikind)%efield)
243 : END IF
244 :
245 14767 : IF (ASSOCIATED(qs_force(ikind)%eev)) THEN
246 9716 : DEALLOCATE (qs_force(ikind)%eev)
247 : END IF
248 : END DO
249 :
250 5051 : DEALLOCATE (qs_force)
251 :
252 5051 : END SUBROUTINE deallocate_qs_force
253 :
254 : ! **************************************************************************************************
255 : !> \brief Initialize a Quickstep force data structure.
256 : !> \param qs_force ...
257 : !> \date 15.07.2002
258 : !> \author MK
259 : !> \version 1.0
260 : ! **************************************************************************************************
261 13663 : SUBROUTINE zero_qs_force(qs_force)
262 :
263 : TYPE(qs_force_type), DIMENSION(:), POINTER :: qs_force
264 :
265 : INTEGER :: ikind
266 :
267 13663 : CPASSERT(ASSOCIATED(qs_force))
268 :
269 39933 : DO ikind = 1, SIZE(qs_force)
270 336866 : qs_force(ikind)%all_potential(:, :) = 0.0_dp
271 336866 : qs_force(ikind)%cneo_potential(:, :) = 0.0_dp
272 336866 : qs_force(ikind)%core_overlap(:, :) = 0.0_dp
273 336866 : qs_force(ikind)%gth_ppl(:, :) = 0.0_dp
274 336866 : qs_force(ikind)%gth_nlcc(:, :) = 0.0_dp
275 336866 : qs_force(ikind)%gth_ppnl(:, :) = 0.0_dp
276 336866 : qs_force(ikind)%kinetic(:, :) = 0.0_dp
277 336866 : qs_force(ikind)%overlap(:, :) = 0.0_dp
278 336866 : qs_force(ikind)%overlap_admm(:, :) = 0.0_dp
279 336866 : qs_force(ikind)%rho_core(:, :) = 0.0_dp
280 336866 : qs_force(ikind)%rho_elec(:, :) = 0.0_dp
281 336866 : qs_force(ikind)%rho_lri_elec(:, :) = 0.0_dp
282 336866 : qs_force(ikind)%rho_cneo_nuc(:, :) = 0.0_dp
283 336866 : qs_force(ikind)%vhxc_atom(:, :) = 0.0_dp
284 336866 : qs_force(ikind)%g0s_Vh_elec(:, :) = 0.0_dp
285 336866 : qs_force(ikind)%repulsive(:, :) = 0.0_dp
286 336866 : qs_force(ikind)%dispersion(:, :) = 0.0_dp
287 336866 : qs_force(ikind)%gcp(:, :) = 0.0_dp
288 336866 : qs_force(ikind)%other(:, :) = 0.0_dp
289 336866 : qs_force(ikind)%fock_4c(:, :) = 0.0_dp
290 336866 : qs_force(ikind)%ehrenfest(:, :) = 0.0_dp
291 336866 : qs_force(ikind)%efield(:, :) = 0.0_dp
292 336866 : qs_force(ikind)%eev(:, :) = 0.0_dp
293 336866 : qs_force(ikind)%mp2_non_sep(:, :) = 0.0_dp
294 350529 : qs_force(ikind)%total(:, :) = 0.0_dp
295 : END DO
296 :
297 13663 : END SUBROUTINE zero_qs_force
298 :
299 : ! **************************************************************************************************
300 : !> \brief Sum up two qs_force entities qs_force_out = qs_force_out + qs_force_in
301 : !> \param qs_force_out ...
302 : !> \param qs_force_in ...
303 : !> \author JGH
304 : ! **************************************************************************************************
305 1588 : SUBROUTINE sum_qs_force(qs_force_out, qs_force_in)
306 :
307 : TYPE(qs_force_type), DIMENSION(:), POINTER :: qs_force_out, qs_force_in
308 :
309 : INTEGER :: ikind
310 :
311 1588 : CPASSERT(ASSOCIATED(qs_force_out))
312 1588 : CPASSERT(ASSOCIATED(qs_force_in))
313 :
314 4856 : DO ikind = 1, SIZE(qs_force_out)
315 : qs_force_out(ikind)%all_potential(:, :) = qs_force_out(ikind)%all_potential(:, :) + &
316 48568 : qs_force_in(ikind)%all_potential(:, :)
317 : qs_force_out(ikind)%cneo_potential(:, :) = qs_force_out(ikind)%cneo_potential(:, :) + &
318 48568 : qs_force_in(ikind)%cneo_potential(:, :)
319 : qs_force_out(ikind)%core_overlap(:, :) = qs_force_out(ikind)%core_overlap(:, :) + &
320 48568 : qs_force_in(ikind)%core_overlap(:, :)
321 : qs_force_out(ikind)%gth_ppl(:, :) = qs_force_out(ikind)%gth_ppl(:, :) + &
322 48568 : qs_force_in(ikind)%gth_ppl(:, :)
323 : qs_force_out(ikind)%gth_nlcc(:, :) = qs_force_out(ikind)%gth_nlcc(:, :) + &
324 48568 : qs_force_in(ikind)%gth_nlcc(:, :)
325 : qs_force_out(ikind)%gth_ppnl(:, :) = qs_force_out(ikind)%gth_ppnl(:, :) + &
326 48568 : qs_force_in(ikind)%gth_ppnl(:, :)
327 : qs_force_out(ikind)%kinetic(:, :) = qs_force_out(ikind)%kinetic(:, :) + &
328 48568 : qs_force_in(ikind)%kinetic(:, :)
329 : qs_force_out(ikind)%overlap(:, :) = qs_force_out(ikind)%overlap(:, :) + &
330 48568 : qs_force_in(ikind)%overlap(:, :)
331 : qs_force_out(ikind)%overlap_admm(:, :) = qs_force_out(ikind)%overlap_admm(:, :) + &
332 48568 : qs_force_in(ikind)%overlap_admm(:, :)
333 : qs_force_out(ikind)%rho_core(:, :) = qs_force_out(ikind)%rho_core(:, :) + &
334 48568 : qs_force_in(ikind)%rho_core(:, :)
335 : qs_force_out(ikind)%rho_elec(:, :) = qs_force_out(ikind)%rho_elec(:, :) + &
336 48568 : qs_force_in(ikind)%rho_elec(:, :)
337 : qs_force_out(ikind)%rho_lri_elec(:, :) = qs_force_out(ikind)%rho_lri_elec(:, :) + &
338 48568 : qs_force_in(ikind)%rho_lri_elec(:, :)
339 : qs_force_out(ikind)%rho_cneo_nuc(:, :) = qs_force_out(ikind)%rho_cneo_nuc(:, :) + &
340 48568 : qs_force_in(ikind)%rho_cneo_nuc(:, :)
341 : qs_force_out(ikind)%vhxc_atom(:, :) = qs_force_out(ikind)%vhxc_atom(:, :) + &
342 48568 : qs_force_in(ikind)%vhxc_atom(:, :)
343 : qs_force_out(ikind)%g0s_Vh_elec(:, :) = qs_force_out(ikind)%g0s_Vh_elec(:, :) + &
344 48568 : qs_force_in(ikind)%g0s_Vh_elec(:, :)
345 : qs_force_out(ikind)%repulsive(:, :) = qs_force_out(ikind)%repulsive(:, :) + &
346 48568 : qs_force_in(ikind)%repulsive(:, :)
347 : qs_force_out(ikind)%dispersion(:, :) = qs_force_out(ikind)%dispersion(:, :) + &
348 48568 : qs_force_in(ikind)%dispersion(:, :)
349 : qs_force_out(ikind)%gcp(:, :) = qs_force_out(ikind)%gcp(:, :) + &
350 48568 : qs_force_in(ikind)%gcp(:, :)
351 : qs_force_out(ikind)%other(:, :) = qs_force_out(ikind)%other(:, :) + &
352 48568 : qs_force_in(ikind)%other(:, :)
353 : qs_force_out(ikind)%fock_4c(:, :) = qs_force_out(ikind)%fock_4c(:, :) + &
354 48568 : qs_force_in(ikind)%fock_4c(:, :)
355 : qs_force_out(ikind)%ehrenfest(:, :) = qs_force_out(ikind)%ehrenfest(:, :) + &
356 48568 : qs_force_in(ikind)%ehrenfest(:, :)
357 : qs_force_out(ikind)%efield(:, :) = qs_force_out(ikind)%efield(:, :) + &
358 48568 : qs_force_in(ikind)%efield(:, :)
359 : qs_force_out(ikind)%eev(:, :) = qs_force_out(ikind)%eev(:, :) + &
360 48568 : qs_force_in(ikind)%eev(:, :)
361 : qs_force_out(ikind)%mp2_non_sep(:, :) = qs_force_out(ikind)%mp2_non_sep(:, :) + &
362 48568 : qs_force_in(ikind)%mp2_non_sep(:, :)
363 : qs_force_out(ikind)%total(:, :) = qs_force_out(ikind)%total(:, :) + &
364 50156 : qs_force_in(ikind)%total(:, :)
365 : END DO
366 :
367 1588 : END SUBROUTINE sum_qs_force
368 :
369 : ! **************************************************************************************************
370 : !> \brief Replicate and sum up the force
371 : !> \param qs_force ...
372 : !> \param para_env ...
373 : !> \date 25.05.2016
374 : !> \author JHU
375 : !> \version 1.0
376 : ! **************************************************************************************************
377 11313 : SUBROUTINE replicate_qs_force(qs_force, para_env)
378 :
379 : TYPE(qs_force_type), DIMENSION(:), POINTER :: qs_force
380 : TYPE(mp_para_env_type), POINTER :: para_env
381 :
382 : INTEGER :: ikind
383 :
384 : ! *** replicate forces ***
385 33231 : DO ikind = 1, SIZE(qs_force)
386 585398 : CALL para_env%sum(qs_force(ikind)%overlap)
387 585398 : CALL para_env%sum(qs_force(ikind)%overlap_admm)
388 585398 : CALL para_env%sum(qs_force(ikind)%kinetic)
389 585398 : CALL para_env%sum(qs_force(ikind)%gth_ppl)
390 585398 : CALL para_env%sum(qs_force(ikind)%gth_nlcc)
391 585398 : CALL para_env%sum(qs_force(ikind)%gth_ppnl)
392 585398 : CALL para_env%sum(qs_force(ikind)%all_potential)
393 585398 : CALL para_env%sum(qs_force(ikind)%cneo_potential)
394 585398 : CALL para_env%sum(qs_force(ikind)%core_overlap)
395 585398 : CALL para_env%sum(qs_force(ikind)%rho_core)
396 585398 : CALL para_env%sum(qs_force(ikind)%rho_elec)
397 585398 : CALL para_env%sum(qs_force(ikind)%rho_lri_elec)
398 585398 : CALL para_env%sum(qs_force(ikind)%rho_cneo_nuc)
399 585398 : CALL para_env%sum(qs_force(ikind)%vhxc_atom)
400 585398 : CALL para_env%sum(qs_force(ikind)%g0s_Vh_elec)
401 585398 : CALL para_env%sum(qs_force(ikind)%fock_4c)
402 585398 : CALL para_env%sum(qs_force(ikind)%mp2_non_sep)
403 585398 : CALL para_env%sum(qs_force(ikind)%repulsive)
404 585398 : CALL para_env%sum(qs_force(ikind)%dispersion)
405 585398 : CALL para_env%sum(qs_force(ikind)%gcp)
406 585398 : CALL para_env%sum(qs_force(ikind)%ehrenfest)
407 :
408 : qs_force(ikind)%total(:, :) = qs_force(ikind)%total(:, :) + &
409 : qs_force(ikind)%core_overlap(:, :) + &
410 : qs_force(ikind)%gth_ppl(:, :) + &
411 : qs_force(ikind)%gth_nlcc(:, :) + &
412 : qs_force(ikind)%gth_ppnl(:, :) + &
413 : qs_force(ikind)%all_potential(:, :) + &
414 : qs_force(ikind)%cneo_potential(:, :) + &
415 : qs_force(ikind)%kinetic(:, :) + &
416 : qs_force(ikind)%overlap(:, :) + &
417 : qs_force(ikind)%overlap_admm(:, :) + &
418 : qs_force(ikind)%rho_core(:, :) + &
419 : qs_force(ikind)%rho_elec(:, :) + &
420 : qs_force(ikind)%rho_lri_elec(:, :) + &
421 : qs_force(ikind)%rho_cneo_nuc(:, :) + &
422 : qs_force(ikind)%vhxc_atom(:, :) + &
423 : qs_force(ikind)%g0s_Vh_elec(:, :) + &
424 : qs_force(ikind)%fock_4c(:, :) + &
425 : qs_force(ikind)%mp2_non_sep(:, :) + &
426 : qs_force(ikind)%repulsive(:, :) + &
427 : qs_force(ikind)%dispersion(:, :) + &
428 : qs_force(ikind)%gcp(:, :) + &
429 : qs_force(ikind)%ehrenfest(:, :) + &
430 : qs_force(ikind)%efield(:, :) + &
431 314971 : qs_force(ikind)%eev(:, :)
432 : END DO
433 :
434 11313 : END SUBROUTINE replicate_qs_force
435 :
436 : ! **************************************************************************************************
437 : !> \brief Add force to a force_type variable.
438 : !> \param force Input force, dimension (3,natom)
439 : !> \param qs_force The force type variable to be used
440 : !> \param forcetype ...
441 : !> \param atomic_kind_set ...
442 : !> \par History
443 : !> 07.2014 JGH
444 : !> \author JGH
445 : ! **************************************************************************************************
446 1248 : SUBROUTINE add_qs_force(force, qs_force, forcetype, atomic_kind_set)
447 :
448 : REAL(KIND=dp), DIMENSION(:, :), INTENT(IN) :: force
449 : TYPE(qs_force_type), DIMENSION(:), POINTER :: qs_force
450 : CHARACTER(LEN=*), INTENT(IN) :: forcetype
451 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
452 :
453 : INTEGER :: ia, iatom, ikind, natom_kind
454 : TYPE(atomic_kind_type), POINTER :: atomic_kind
455 :
456 : ! ------------------------------------------------------------------------
457 :
458 1248 : CPASSERT(ASSOCIATED(qs_force))
459 :
460 1248 : SELECT CASE (forcetype)
461 : CASE ("overlap_admm")
462 3532 : DO ikind = 1, SIZE(atomic_kind_set, 1)
463 2284 : atomic_kind => atomic_kind_set(ikind)
464 2284 : CALL get_atomic_kind(atomic_kind=atomic_kind, natom=natom_kind)
465 7128 : DO ia = 1, natom_kind
466 3596 : iatom = atomic_kind%atom_list(ia)
467 16668 : qs_force(ikind)%overlap_admm(:, ia) = qs_force(ikind)%overlap_admm(:, ia) + force(:, iatom)
468 : END DO
469 : END DO
470 : CASE DEFAULT
471 : CALL cp_abort(__LOCATION__, &
472 : "<overlap_admm> is supported as the <forcetype> "// &
473 : "for add_qs_force, found unknown option "// &
474 1248 : "<"//TRIM(forcetype)//">")
475 : END SELECT
476 :
477 1248 : END SUBROUTINE add_qs_force
478 :
479 : ! **************************************************************************************************
480 : !> \brief Put force to a force_type variable.
481 : !> \param force Input force, dimension (3,natom)
482 : !> \param qs_force The force type variable to be used
483 : !> \param forcetype ...
484 : !> \param atomic_kind_set ...
485 : !> \par History
486 : !> 09.2019 JGH
487 : !> \author JGH
488 : ! **************************************************************************************************
489 0 : SUBROUTINE put_qs_force(force, qs_force, forcetype, atomic_kind_set)
490 :
491 : REAL(KIND=dp), DIMENSION(:, :), INTENT(IN) :: force
492 : TYPE(qs_force_type), DIMENSION(:), POINTER :: qs_force
493 : CHARACTER(LEN=*), INTENT(IN) :: forcetype
494 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
495 :
496 : INTEGER :: ia, iatom, ikind, natom_kind
497 : TYPE(atomic_kind_type), POINTER :: atomic_kind
498 :
499 : ! ------------------------------------------------------------------------
500 :
501 0 : SELECT CASE (forcetype)
502 : CASE ("dispersion")
503 0 : DO ikind = 1, SIZE(atomic_kind_set, 1)
504 0 : atomic_kind => atomic_kind_set(ikind)
505 0 : CALL get_atomic_kind(atomic_kind=atomic_kind, natom=natom_kind)
506 0 : DO ia = 1, natom_kind
507 0 : iatom = atomic_kind%atom_list(ia)
508 0 : qs_force(ikind)%dispersion(:, ia) = force(:, iatom)
509 : END DO
510 : END DO
511 : CASE DEFAULT
512 : CALL cp_abort(__LOCATION__, &
513 : "<dispersion> is supported as the <forcetype> "// &
514 : "for put_qs_force, found unknown option "// &
515 0 : "<"//TRIM(forcetype)//">")
516 : END SELECT
517 :
518 0 : END SUBROUTINE put_qs_force
519 :
520 : ! **************************************************************************************************
521 : !> \brief Get force from a force_type variable.
522 : !> \param force Input force, dimension (3,natom)
523 : !> \param qs_force The force type variable to be used
524 : !> \param forcetype ...
525 : !> \param atomic_kind_set ...
526 : !> \par History
527 : !> 09.2019 JGH
528 : !> \author JGH
529 : ! **************************************************************************************************
530 0 : SUBROUTINE get_qs_force(force, qs_force, forcetype, atomic_kind_set)
531 :
532 : REAL(KIND=dp), DIMENSION(:, :), INTENT(INOUT) :: force
533 : TYPE(qs_force_type), DIMENSION(:), POINTER :: qs_force
534 : CHARACTER(LEN=*), INTENT(IN) :: forcetype
535 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
536 :
537 : INTEGER :: ia, iatom, ikind, natom_kind
538 : TYPE(atomic_kind_type), POINTER :: atomic_kind
539 :
540 : ! ------------------------------------------------------------------------
541 :
542 0 : SELECT CASE (forcetype)
543 : CASE ("dispersion")
544 0 : DO ikind = 1, SIZE(atomic_kind_set, 1)
545 0 : atomic_kind => atomic_kind_set(ikind)
546 0 : CALL get_atomic_kind(atomic_kind=atomic_kind, natom=natom_kind)
547 0 : DO ia = 1, natom_kind
548 0 : iatom = atomic_kind%atom_list(ia)
549 0 : force(:, iatom) = qs_force(ikind)%dispersion(:, ia)
550 : END DO
551 : END DO
552 : CASE DEFAULT
553 : CALL cp_abort(__LOCATION__, &
554 : "<dispersion> is supported as the <forcetype> "// &
555 : "for get_qs_force, found unknown option "// &
556 0 : "<"//TRIM(forcetype)//">")
557 : END SELECT
558 :
559 0 : END SUBROUTINE get_qs_force
560 :
561 : ! **************************************************************************************************
562 : !> \brief Get current total force
563 : !> \param force Input force, dimension (3,natom)
564 : !> \param qs_force The force type variable to be used
565 : !> \param atomic_kind_set ...
566 : !> \par History
567 : !> 09.2019 JGH
568 : !> \author JGH
569 : ! **************************************************************************************************
570 998 : SUBROUTINE total_qs_force(force, qs_force, atomic_kind_set)
571 :
572 : REAL(KIND=dp), DIMENSION(:, :), INTENT(INOUT) :: force
573 : TYPE(qs_force_type), DIMENSION(:), POINTER :: qs_force
574 : TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
575 :
576 : INTEGER :: ia, iatom, ikind, natom_kind
577 : TYPE(atomic_kind_type), POINTER :: atomic_kind
578 :
579 : ! ------------------------------------------------------------------------
580 :
581 13230 : force(:, :) = 0.0_dp
582 3062 : DO ikind = 1, SIZE(atomic_kind_set, 1)
583 2064 : atomic_kind => atomic_kind_set(ikind)
584 2064 : CALL get_atomic_kind(atomic_kind=atomic_kind, natom=natom_kind)
585 6120 : DO ia = 1, natom_kind
586 3058 : iatom = atomic_kind%atom_list(ia)
587 : force(:, iatom) = qs_force(ikind)%core_overlap(:, ia) + &
588 : qs_force(ikind)%gth_ppl(:, ia) + &
589 : qs_force(ikind)%gth_nlcc(:, ia) + &
590 : qs_force(ikind)%gth_ppnl(:, ia) + &
591 : qs_force(ikind)%all_potential(:, ia) + &
592 : qs_force(ikind)%cneo_potential(:, ia) + &
593 : qs_force(ikind)%kinetic(:, ia) + &
594 : qs_force(ikind)%overlap(:, ia) + &
595 : qs_force(ikind)%overlap_admm(:, ia) + &
596 : qs_force(ikind)%rho_core(:, ia) + &
597 : qs_force(ikind)%rho_elec(:, ia) + &
598 : qs_force(ikind)%rho_lri_elec(:, ia) + &
599 : qs_force(ikind)%rho_cneo_nuc(:, ia) + &
600 : qs_force(ikind)%vhxc_atom(:, ia) + &
601 : qs_force(ikind)%g0s_Vh_elec(:, ia) + &
602 : qs_force(ikind)%fock_4c(:, ia) + &
603 : qs_force(ikind)%mp2_non_sep(:, ia) + &
604 : qs_force(ikind)%repulsive(:, ia) + &
605 : qs_force(ikind)%dispersion(:, ia) + &
606 : qs_force(ikind)%gcp(:, ia) + &
607 : qs_force(ikind)%ehrenfest(:, ia) + &
608 : qs_force(ikind)%efield(:, ia) + &
609 14296 : qs_force(ikind)%eev(:, ia)
610 : END DO
611 : END DO
612 :
613 998 : END SUBROUTINE total_qs_force
614 :
615 : ! **************************************************************************************************
616 : !> \brief Write a Quickstep force data for 1 atom
617 : !> \param qs_force ...
618 : !> \param ikind ...
619 : !> \param iatom ...
620 : !> \param iunit ...
621 : !> \date 05.06.2002
622 : !> \author MK/JGH
623 : !> \version 1.0
624 : ! **************************************************************************************************
625 0 : SUBROUTINE write_forces_debug(qs_force, ikind, iatom, iunit)
626 :
627 : TYPE(qs_force_type), DIMENSION(:), POINTER :: qs_force
628 : INTEGER, INTENT(IN), OPTIONAL :: ikind, iatom, iunit
629 :
630 : CHARACTER(LEN=35) :: fmtstr2
631 : CHARACTER(LEN=48) :: fmtstr1
632 : INTEGER :: iounit, jatom, jkind
633 : REAL(KIND=dp), DIMENSION(3) :: total
634 : TYPE(cp_logger_type), POINTER :: logger
635 :
636 0 : IF (PRESENT(iunit)) THEN
637 0 : iounit = iunit
638 : ELSE
639 0 : NULLIFY (logger)
640 0 : logger => cp_get_default_logger()
641 0 : iounit = cp_logger_get_default_io_unit(logger)
642 : END IF
643 0 : IF (PRESENT(ikind)) THEN
644 0 : jkind = ikind
645 : ELSE
646 0 : jkind = 1
647 : END IF
648 0 : IF (PRESENT(iatom)) THEN
649 0 : jatom = iatom
650 : ELSE
651 0 : jatom = 1
652 : END IF
653 :
654 0 : IF (iounit > 0) THEN
655 :
656 0 : fmtstr1 = "(/,T2,A,/,T3,A,T11,A,T23,A,T40,A1,2(17X,A1))"
657 0 : fmtstr2 = "((T2,I5,4X,I4,T18,A,T34,3F18.12))"
658 :
659 : WRITE (UNIT=iounit, FMT=fmtstr1) &
660 0 : "FORCES [a.u.]", "Atom", "Kind", "Component", "X", "Y", "Z"
661 :
662 : total(1:3) = qs_force(jkind)%overlap(1:3, jatom) &
663 : + qs_force(jkind)%overlap_admm(1:3, jatom) &
664 : + qs_force(jkind)%kinetic(1:3, jatom) &
665 : + qs_force(jkind)%gth_ppl(1:3, jatom) &
666 : + qs_force(jkind)%gth_ppnl(1:3, jatom) &
667 : + qs_force(jkind)%gth_nlcc(1:3, jatom) &
668 : + qs_force(jkind)%all_potential(1:3, jatom) &
669 : + qs_force(jkind)%cneo_potential(1:3, jatom) &
670 : + qs_force(jkind)%rho_cneo_nuc(1:3, jatom) &
671 : + qs_force(jkind)%core_overlap(1:3, jatom) &
672 : + qs_force(jkind)%rho_core(1:3, jatom) &
673 : + qs_force(jkind)%rho_elec(1:3, jatom) &
674 : + qs_force(jkind)%rho_lri_elec(1:3, jatom) &
675 : + qs_force(jkind)%vhxc_atom(1:3, jatom) &
676 : + qs_force(jkind)%g0s_Vh_elec(1:3, jatom) &
677 : + qs_force(jkind)%dispersion(1:3, jatom) &
678 : + qs_force(jkind)%repulsive(1:3, jatom) &
679 : + qs_force(jkind)%gcp(1:3, jatom) &
680 : + qs_force(jkind)%efield(1:3, jatom) &
681 : + qs_force(jkind)%eev(1:3, jatom) &
682 : + qs_force(jkind)%ehrenfest(1:3, jatom) &
683 : + qs_force(jkind)%fock_4c(1:3, jatom) &
684 0 : + qs_force(jkind)%mp2_non_sep(1:3, jatom)
685 :
686 : WRITE (UNIT=iounit, FMT=fmtstr2) &
687 0 : jatom, jkind, " overlap", qs_force(jkind)%overlap(1:3, jatom), &
688 0 : jatom, jkind, " overlap_admm", qs_force(jkind)%overlap_admm(1:3, jatom), &
689 0 : jatom, jkind, " kinetic", qs_force(jkind)%kinetic(1:3, jatom), &
690 0 : jatom, jkind, " gth_ppl", qs_force(jkind)%gth_ppl(1:3, jatom), &
691 0 : jatom, jkind, " gth_ppnl", qs_force(jkind)%gth_ppnl(1:3, jatom), &
692 0 : jatom, jkind, " gth_nlcc", qs_force(jkind)%gth_nlcc(1:3, jatom), &
693 0 : jatom, jkind, " all_potential", qs_force(jkind)%all_potential(1:3, jatom), &
694 0 : jatom, jkind, "cneo_potential", qs_force(jkind)%cneo_potential(1:3, jatom), &
695 0 : jatom, jkind, " rho_cneo_nuc", qs_force(jkind)%rho_cneo_nuc(1:3, jatom), &
696 0 : jatom, jkind, " core_overlap", qs_force(jkind)%core_overlap(1:3, jatom), &
697 0 : jatom, jkind, " rho_core", qs_force(jkind)%rho_core(1:3, jatom), &
698 0 : jatom, jkind, " rho_elec", qs_force(jkind)%rho_elec(1:3, jatom), &
699 0 : jatom, jkind, " rho_lri_elec", qs_force(jkind)%rho_lri_elec(1:3, jatom), &
700 0 : jatom, jkind, " vhxc_atom", qs_force(jkind)%vhxc_atom(1:3, jatom), &
701 0 : jatom, jkind, " g0s_Vh_elec", qs_force(jkind)%g0s_Vh_elec(1:3, jatom), &
702 0 : jatom, jkind, " dispersion", qs_force(jkind)%dispersion(1:3, jatom), &
703 0 : jatom, jkind, " repulsive", qs_force(jkind)%repulsive(1:3, jatom), &
704 0 : jatom, jkind, " gcp", qs_force(jkind)%gcp(1:3, jatom), &
705 0 : jatom, jkind, " efield", qs_force(jkind)%efield(1:3, jatom), &
706 0 : jatom, jkind, " eev", qs_force(jkind)%eev(1:3, jatom), &
707 0 : jatom, jkind, " ehrenfest", qs_force(jkind)%ehrenfest(1:3, jatom), &
708 0 : jatom, jkind, " fock_4c", qs_force(jkind)%fock_4c(1:3, jatom), &
709 0 : jatom, jkind, " mp2_non_sep", qs_force(jkind)%mp2_non_sep(1:3, jatom), &
710 0 : jatom, jkind, " total", total(1:3)
711 :
712 : END IF
713 :
714 0 : END SUBROUTINE write_forces_debug
715 :
716 0 : END MODULE qs_force_types
|
