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