(git:374b731)
Loading...
Searching...
No Matches
fist_force.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!> Torsions added(DG)05-Dec-2000
11!> Variable names changed(DG)05-Dec-2000
12!> CJM SEPT-12-2002: int_env is now passed
13!> CJM NOV-30-2003: only uses fist_env
14!> \author CJM & JGH
15! **************************************************************************************************
19 USE atprop_types, ONLY: atprop_type
20 USE cell_methods, ONLY: init_cell
21 USE cell_types, ONLY: cell_type
24 USE cp_output_handling, ONLY: cp_p_file,&
36 USE ewalds, ONLY: ewald_evaluate,&
57 USE kinds, ONLY: dp
66 USE pme, ONLY: pme_evaluate
72 USE spme, ONLY: spme_evaluate
73 USE virial_types, ONLY: virial_type
74#include "./base/base_uses.f90"
75
76 IMPLICIT NONE
77 PRIVATE
79 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'fist_force'
80
81 TYPE debug_variables_type
82 REAL(KIND=dp) :: pot_manybody, pot_bend, pot_torsion
83 REAL(KIND=dp) :: pot_nonbond, pot_g, pot_bond
84 REAL(KIND=dp) :: pot_imptors, pot_urey_bradley
85 REAL(KIND=dp) :: pot_opbend
86 REAL(KIND=dp), DIMENSION(:, :), POINTER :: f_nonbond, f_g, f_bond, f_bend, &
87 f_torsion, f_imptors, f_ub, &
88 f_opbend
89 REAL(KIND=dp), DIMENSION(3, 3) :: pv_nonbond, pv_g, pv_bond, pv_ub, &
90 pv_bend, pv_torsion, pv_imptors, &
91 pv_opbend
92 END TYPE debug_variables_type
93
94CONTAINS
95
96! **************************************************************************************************
97!> \brief Calculates the total potential energy, total force, and the
98!> total pressure tensor from the potentials
99!> \param fist_env ...
100!> \param debug ...
101!> \par History
102!> Harald Forbert(Dec-2000): Changes for multiple linked lists
103!> cjm, 20-Feb-2001: box_ref used to initialize ewald. Now
104!> have consistent restarts with npt and ewald
105!> JGH(15-Mar-2001): box_change replaces ensemble parameter
106!> Call ewald_setup if box has changed
107!> Consistent setup for PME and SPME
108!> cjm, 28-Feb-2006: box_change is gone
109!> \author CJM & JGH
110! **************************************************************************************************
111 SUBROUTINE fist_calc_energy_force(fist_env, debug)
112 TYPE(fist_environment_type), POINTER :: fist_env
113 TYPE(debug_variables_type), INTENT(INOUT), &
114 OPTIONAL :: debug
115
116 CHARACTER(len=*), PARAMETER :: routinen = 'fist_calc_energy_force'
117
118 INTEGER :: do_ipol, ewald_type, fg_coulomb_size, handle, i, ii, ikind, iparticle_kind, &
119 iparticle_local, iw, iw2, j, natoms, nlocal_particles, node, nparticle_kind, &
120 nparticle_local, nshell, shell_index
121 LOGICAL :: do_multipoles, shell_model_ad, &
122 shell_present, use_virial
123 REAL(kind=dp) :: ef_ener, fc, fs, mass, pot_bend, pot_bond, pot_imptors, pot_manybody, &
124 pot_nonbond, pot_opbend, pot_shell, pot_torsion, pot_urey_bradley, vg_coulomb, xdum1, &
125 xdum2, xdum3
126 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: ef_f, f_nonbond, f_total, fcore_nonbond, &
127 fcore_shell_bonded, fcore_total, fg_coulomb, fgcore_coulomb, fgshell_coulomb, &
128 fshell_nonbond, fshell_total
129 REAL(kind=dp), DIMENSION(3, 3) :: ef_pv, ident, pv_bc, pv_bend, pv_bond, &
130 pv_g, pv_imptors, pv_nonbond, &
131 pv_opbend, pv_torsion, pv_urey_bradley
132 REAL(kind=dp), DIMENSION(:), POINTER :: e_coulomb
133 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: pv_coulomb
134 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
135 TYPE(atomic_kind_type), POINTER :: atomic_kind
136 TYPE(atprop_type), POINTER :: atprop_env
137 TYPE(cell_type), POINTER :: cell
138 TYPE(cp_logger_type), POINTER :: logger
139 TYPE(cp_subsys_type), POINTER :: subsys
140 TYPE(distribution_1d_type), POINTER :: local_molecules, local_particles
141 TYPE(ewald_environment_type), POINTER :: ewald_env
142 TYPE(ewald_pw_type), POINTER :: ewald_pw
143 TYPE(exclusion_type), DIMENSION(:), POINTER :: exclusions
144 TYPE(fist_efield_type), POINTER :: efield
145 TYPE(fist_energy_type), POINTER :: thermo
146 TYPE(fist_nonbond_env_type), POINTER :: fist_nonbond_env
147 TYPE(molecule_kind_type), DIMENSION(:), POINTER :: molecule_kind_set
148 TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set
149 TYPE(mp_para_env_type), POINTER :: para_env
150 TYPE(multipole_type), POINTER :: multipoles
151 TYPE(particle_type), DIMENSION(:), POINTER :: core_particle_set, particle_set, &
152 shell_particle_set
153 TYPE(section_vals_type), POINTER :: force_env_section, mm_section, &
154 print_section
155 TYPE(shell_kind_type), POINTER :: shell
156 TYPE(virial_type), POINTER :: virial
157
158 CALL timeset(routinen, handle)
159 NULLIFY (logger)
160 NULLIFY (subsys, virial, atprop_env, para_env, force_env_section)
161 logger => cp_get_default_logger()
162
163 CALL fist_env_get(fist_env, &
164 subsys=subsys, &
165 para_env=para_env, &
166 input=force_env_section)
167
168 CALL cp_subsys_get(subsys, &
169 virial=virial, &
170 atprop=atprop_env)
171
172 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
173 NULLIFY (atomic_kind, atomic_kind_set, cell, ewald_pw, ewald_env, &
174 fist_nonbond_env, mm_section, local_molecules, local_particles, &
175 molecule_kind_set, molecule_set, particle_set, print_section, &
176 shell, shell_particle_set, core_particle_set, thermo, multipoles, &
177 e_coulomb, pv_coulomb)
178
179 mm_section => section_vals_get_subs_vals(force_env_section, "MM")
180 iw = cp_print_key_unit_nr(logger, mm_section, "PRINT%DERIVATIVES", &
181 extension=".mmLog")
182 iw2 = cp_print_key_unit_nr(logger, mm_section, "PRINT%EWALD_INFO", &
183 extension=".mmLog")
184
185 CALL fist_env_get(fist_env, ewald_pw=ewald_pw, ewald_env=ewald_env, &
186 local_particles=local_particles, particle_set=particle_set, &
187 atomic_kind_set=atomic_kind_set, molecule_set=molecule_set, &
188 local_molecules=local_molecules, thermo=thermo, &
189 molecule_kind_set=molecule_kind_set, fist_nonbond_env=fist_nonbond_env, &
190 cell=cell, shell_model=shell_present, shell_model_ad=shell_model_ad, &
191 multipoles=multipoles, exclusions=exclusions, efield=efield)
192
193 CALL ewald_env_get(ewald_env, ewald_type=ewald_type, do_multipoles=do_multipoles, &
194 do_ipol=do_ipol)
195
196 ! Initialize ewald grids
197 CALL init_cell(cell)
198 CALL ewald_pw_grid_update(ewald_pw, ewald_env, cell%hmat)
199
200 natoms = SIZE(particle_set)
201 nlocal_particles = 0
202 nparticle_kind = SIZE(atomic_kind_set)
203 DO ikind = 1, nparticle_kind
204 nlocal_particles = nlocal_particles + local_particles%n_el(ikind)
205 END DO
206
207 ALLOCATE (f_nonbond(3, natoms))
208 f_nonbond = 0.0_dp
209
210 nshell = 0
211 IF (shell_present) THEN
212 CALL fist_env_get(fist_env, shell_particle_set=shell_particle_set, &
213 core_particle_set=core_particle_set)
214 cpassert(ASSOCIATED(shell_particle_set))
215 cpassert(ASSOCIATED(core_particle_set))
216 nshell = SIZE(shell_particle_set)
217 ALLOCATE (fshell_nonbond(3, nshell))
218 fshell_nonbond = 0.0_dp
219 ALLOCATE (fcore_nonbond(3, nshell))
220 fcore_nonbond = 0.0_dp
221 ELSE
222 NULLIFY (shell_particle_set, core_particle_set)
223 END IF
224
225 IF (fist_nonbond_env%do_nonbonded) THEN
226 ! First check with list_control to update neighbor lists
227 IF (ASSOCIATED(exclusions)) THEN
228 CALL list_control(atomic_kind_set, particle_set, local_particles, &
229 cell, fist_nonbond_env, para_env, mm_section, shell_particle_set, &
230 core_particle_set, exclusions=exclusions)
231 ELSE
232 CALL list_control(atomic_kind_set, particle_set, local_particles, &
233 cell, fist_nonbond_env, para_env, mm_section, shell_particle_set, &
234 core_particle_set)
235 END IF
236 END IF
237
238 ! Initialize force, energy and pressure tensor arrays
239 DO i = 1, natoms
240 particle_set(i)%f(1) = 0.0_dp
241 particle_set(i)%f(2) = 0.0_dp
242 particle_set(i)%f(3) = 0.0_dp
243 END DO
244 IF (nshell > 0) THEN
245 DO i = 1, nshell
246 shell_particle_set(i)%f(1) = 0.0_dp
247 shell_particle_set(i)%f(2) = 0.0_dp
248 shell_particle_set(i)%f(3) = 0.0_dp
249 core_particle_set(i)%f(1) = 0.0_dp
250 core_particle_set(i)%f(2) = 0.0_dp
251 core_particle_set(i)%f(3) = 0.0_dp
252 END DO
253 END IF
254
255 IF (use_virial) THEN
256 pv_bc = 0.0_dp
257 pv_bond = 0.0_dp
258 pv_bend = 0.0_dp
259 pv_torsion = 0.0_dp
260 pv_imptors = 0.0_dp
261 pv_opbend = 0.0_dp
262 pv_urey_bradley = 0.0_dp
263 pv_nonbond = 0.0_dp
264 pv_g = 0.0_dp
265 virial%pv_virial = 0.0_dp
266 END IF
267
268 pot_nonbond = 0.0_dp
269 pot_manybody = 0.0_dp
270 pot_bond = 0.0_dp
271 pot_bend = 0.0_dp
272 pot_torsion = 0.0_dp
273 pot_opbend = 0.0_dp
274 pot_imptors = 0.0_dp
275 pot_urey_bradley = 0.0_dp
276 pot_shell = 0.0_dp
277 vg_coulomb = 0.0_dp
278 thermo%pot = 0.0_dp
279 thermo%harm_shell = 0.0_dp
280
281 ! Get real-space non-bonded forces:
282 IF (iw > 0) THEN
283 WRITE (iw, '(A)') " FIST:: FORCES IN INPUT..."
284 WRITE (iw, '(3f15.9)') ((particle_set(i)%f(j), j=1, 3), i=1, SIZE(particle_set))
285 END IF
286
287 IF (fist_nonbond_env%do_nonbonded) THEN
288 ! Compute density for EAM
289 CALL density_nonbond(fist_nonbond_env, particle_set, cell, para_env)
290
291 ! Compute embedding function and manybody energy
292 CALL energy_manybody(fist_nonbond_env, atomic_kind_set, local_particles, particle_set, &
293 cell, pot_manybody, para_env, mm_section)
294
295 ! Nonbond contribution + Manybody Forces
296 IF (shell_present) THEN
297 CALL force_nonbond(fist_nonbond_env, ewald_env, particle_set, cell, &
298 pot_nonbond, f_nonbond, pv_nonbond, &
299 fshell_nonbond=fshell_nonbond, fcore_nonbond=fcore_nonbond, &
300 atprop_env=atprop_env, &
301 atomic_kind_set=atomic_kind_set, use_virial=use_virial)
302 ELSE
303 CALL force_nonbond(fist_nonbond_env, ewald_env, particle_set, cell, &
304 pot_nonbond, f_nonbond, pv_nonbond, atprop_env=atprop_env, &
305 atomic_kind_set=atomic_kind_set, use_virial=use_virial)
306 CALL force_nonbond_manybody(fist_nonbond_env, particle_set, cell, f_nonbond, pv_nonbond, &
307 use_virial=use_virial)
308 END IF
309 END IF
310
311 IF (iw > 0) THEN
312 WRITE (iw, '(A)') " FIST:: NONBOND + R-SPACE ELECTROSTATIC FORCES ..."
313 WRITE (iw, '(3f15.9)') f_nonbond
314 IF (shell_present .AND. shell_model_ad) THEN
315 WRITE (iw, '(A)') " FIST:: SHELL NONBOND + R-SPACE ELECTROSTATIC FORCES ..."
316 WRITE (iw, '(3f15.9)') fshell_nonbond
317 WRITE (iw, '(A)') " FIST:: CORE NONBOND + R-SPACE ELECTROSTATIC FORCES ..."
318 WRITE (iw, '(3f15.9)') fcore_nonbond
319 END IF
320 END IF
321
322 ! Get g-space non-bonded forces:
323 IF (ewald_type /= do_ewald_none) THEN
324 ! Determine size of the forces array
325 SELECT CASE (ewald_type)
326 CASE (do_ewald_ewald)
327 fg_coulomb_size = nlocal_particles
328 CASE DEFAULT
329 fg_coulomb_size = natoms
330 END SELECT
331 ! Allocate and zeroing arrays
332 ALLOCATE (fg_coulomb(3, fg_coulomb_size))
333 fg_coulomb = 0.0_dp
334 IF (shell_present) THEN
335 ALLOCATE (fgshell_coulomb(3, nshell))
336 ALLOCATE (fgcore_coulomb(3, nshell))
337 fgshell_coulomb = 0.0_dp
338 fgcore_coulomb = 0.0_dp
339 END IF
340 IF (shell_present .AND. do_multipoles) THEN
341 cpabort("Multipoles and Core-Shell model not implemented.")
342 END IF
343 ! If not multipole: Compute self-interaction and neutralizing background
344 ! for multipoles is handled separately..
345 IF (.NOT. do_multipoles) THEN
346 CALL ewald_self(ewald_env, cell, atomic_kind_set, local_particles, &
347 thermo%e_self, thermo%e_neut, fist_nonbond_env%charges)
348 IF (atprop_env%energy) THEN
349 CALL ewald_self_atom(ewald_env, atomic_kind_set, local_particles, &
350 atprop_env%atener, fist_nonbond_env%charges)
351 atprop_env%atener = atprop_env%atener + thermo%e_neut/SIZE(atprop_env%atener)
352 END IF
353 END IF
354
355 ! Polarizable force-field
356 IF (do_ipol /= do_fist_pol_none) THEN
357 ! Check if an array of chagres was provided and in case abort due to lack of implementation
358 IF (ASSOCIATED(fist_nonbond_env%charges)) THEN
359 cpabort("Polarizable force field and array charges not implemented!")
360 END IF
361 ! Converge the dipoles self-consistently
362 CALL fist_pol_evaluate(atomic_kind_set, multipoles, ewald_env, &
363 ewald_pw, fist_nonbond_env, cell, particle_set, &
364 local_particles, thermo, vg_coulomb, pot_nonbond, f_nonbond, &
365 fg_coulomb, use_virial, pv_g, pv_nonbond, mm_section, do_ipol)
366 ELSE
367 ! Non-Polarizable force-field
368 SELECT CASE (ewald_type)
369 CASE (do_ewald_ewald)
370 ! Parallelisation over atoms --> allocate local atoms
371 IF (shell_present) THEN
372 ! Check if an array of chagres was provided and in case abort due to lack of implementation
373 IF (ASSOCIATED(fist_nonbond_env%charges)) THEN
374 cpabort("Core-Shell and array charges not implemented!")
375 END IF
376 IF (do_multipoles) THEN
377 cpabort("Multipole Ewald and CORE-SHELL not yet implemented within Ewald sum!")
378 ELSE
379 cpabort("Core-Shell model not yet implemented within Ewald sums.")
380 END IF
381 ELSE
382 IF (do_multipoles) THEN
383 ! Check if an array of chagres was provided and in case abort due to lack of implementation
384 IF (ASSOCIATED(fist_nonbond_env%charges)) THEN
385 cpabort("Multipole Ewald and array charges not implemented!")
386 END IF
387 CALL ewald_multipole_evaluate(ewald_env, ewald_pw, fist_nonbond_env, cell, &
388 particle_set, local_particles, vg_coulomb, pot_nonbond, thermo%e_neut, &
389 thermo%e_self, multipoles%task, do_correction_bonded=.true., do_forces=.true., &
390 do_stress=use_virial, do_efield=.false., radii=multipoles%radii, &
391 charges=multipoles%charges, dipoles=multipoles%dipoles, &
392 quadrupoles=multipoles%quadrupoles, forces_local=fg_coulomb, &
393 forces_glob=f_nonbond, pv_local=pv_g, pv_glob=pv_nonbond, iw=iw2, &
394 do_debug=.true., atomic_kind_set=atomic_kind_set, mm_section=mm_section)
395 ELSE
396 IF (atprop_env%energy) THEN
397 ALLOCATE (e_coulomb(fg_coulomb_size))
398 END IF
399 IF (atprop_env%stress) THEN
400 ALLOCATE (pv_coulomb(3, 3, fg_coulomb_size))
401 END IF
402 CALL ewald_evaluate(ewald_env, ewald_pw, cell, atomic_kind_set, particle_set, &
403 local_particles, fg_coulomb, vg_coulomb, pv_g, use_virial=use_virial, &
404 charges=fist_nonbond_env%charges, e_coulomb=e_coulomb, &
405 pv_coulomb=pv_coulomb)
406 END IF
407 END IF
408 CASE (do_ewald_pme)
409 ! Parallelisation over grids --> allocate all atoms
410 IF (shell_present) THEN
411 ! Check if an array of chagres was provided and in case abort due to lack of implementation
412 IF (ASSOCIATED(fist_nonbond_env%charges)) THEN
413 cpabort("Core-Shell and array charges not implemented!")
414 END IF
415 IF (do_multipoles) THEN
416 cpabort("Multipole Ewald and CORE-SHELL not yet implemented within a PME scheme!")
417 ELSE
418 CALL pme_evaluate(ewald_env, ewald_pw, cell, particle_set, vg_coulomb, fg_coulomb, &
419 pv_g, shell_particle_set=shell_particle_set, core_particle_set=core_particle_set, &
420 fgshell_coulomb=fgshell_coulomb, fgcore_coulomb=fgcore_coulomb, use_virial=use_virial, &
421 atprop=atprop_env)
422 CALL para_env%sum(fgshell_coulomb)
423 CALL para_env%sum(fgcore_coulomb)
424 END IF
425 ELSE
426 IF (do_multipoles) THEN
427 cpabort("Multipole Ewald not yet implemented within a PME scheme!")
428 ELSE
429 CALL pme_evaluate(ewald_env, ewald_pw, cell, particle_set, vg_coulomb, fg_coulomb, &
430 pv_g, use_virial=use_virial, charges=fist_nonbond_env%charges, &
431 atprop=atprop_env)
432 END IF
433 END IF
434 CALL para_env%sum(fg_coulomb)
435 CASE (do_ewald_spme)
436 ! Parallelisation over grids --> allocate all atoms
437 IF (shell_present) THEN
438 ! Check if an array of charges was provided and in case abort due to lack of implementation
439 IF (ASSOCIATED(fist_nonbond_env%charges)) THEN
440 cpabort("Core-Shell and array charges not implemented!")
441 END IF
442 IF (do_multipoles) THEN
443 cpabort("Multipole Ewald and CORE-SHELL not yet implemented within a SPME scheme!")
444 ELSE
445 CALL spme_evaluate(ewald_env, ewald_pw, cell, particle_set, fg_coulomb, vg_coulomb, &
446 pv_g, shell_particle_set=shell_particle_set, core_particle_set=core_particle_set, &
447 fgshell_coulomb=fgshell_coulomb, fgcore_coulomb=fgcore_coulomb, use_virial=use_virial, &
448 atprop=atprop_env)
449 CALL para_env%sum(fgshell_coulomb)
450 CALL para_env%sum(fgcore_coulomb)
451 END IF
452 ELSE
453 IF (do_multipoles) THEN
454 cpabort("Multipole Ewald not yet implemented within a SPME scheme!")
455 ELSE
456 CALL spme_evaluate(ewald_env, ewald_pw, cell, particle_set, fg_coulomb, vg_coulomb, &
457 pv_g, use_virial=use_virial, charges=fist_nonbond_env%charges, &
458 atprop=atprop_env)
459 END IF
460 END IF
461 CALL para_env%sum(fg_coulomb)
462 END SELECT
463 END IF
464
465 ! Subtract the interaction between screening charges. This is a
466 ! correction in real-space and depends on the neighborlists. Therefore
467 ! it is only executed if fist_nonbond_env%do_nonbonded is set.
468 IF (fist_nonbond_env%do_nonbonded) THEN
469 ! Correction for core-shell model
470 IF (shell_present) THEN
471 CALL bonded_correct_gaussian(fist_nonbond_env, atomic_kind_set, &
472 local_particles, particle_set, ewald_env, thermo%e_bonded, &
473 pv_bc, shell_particle_set=shell_particle_set, &
474 core_particle_set=core_particle_set, atprop_env=atprop_env, cell=cell, &
475 use_virial=use_virial)
476 ELSE
477 IF (.NOT. do_multipoles) THEN
478 CALL bonded_correct_gaussian(fist_nonbond_env, &
479 atomic_kind_set, local_particles, particle_set, &
480 ewald_env, thermo%e_bonded, pv_bc=pv_bc, atprop_env=atprop_env, cell=cell, &
481 use_virial=use_virial)
482 END IF
483 END IF
484 END IF
485
486 IF (.NOT. do_multipoles) THEN
487 ! Multipole code has its own printing routine.
488 CALL ewald_print(iw2, pot_nonbond, vg_coulomb, thermo%e_self, thermo%e_neut, &
489 thermo%e_bonded)
490 END IF
491 ELSE
492 IF (use_virial) THEN
493 pv_g = 0.0_dp
494 pv_bc = 0.0_dp
495 END IF
496 thermo%e_neut = 0.0_dp
497 END IF
498
499 IF (iw > 0) THEN
500 IF (ALLOCATED(fg_coulomb)) THEN
501 WRITE (iw, '(A)') " FIST:: NONBONDED ELECTROSTATIC FORCES IN G-SPACE..."
502 WRITE (iw, '(3f15.9)') ((fg_coulomb(j, i), j=1, 3), i=1, SIZE(fg_coulomb, 2))
503 END IF
504 IF (shell_present .AND. shell_model_ad .AND. ALLOCATED(fgshell_coulomb)) THEN
505 WRITE (iw, '(A)') " FIST:: SHELL NONBONDED ELECTROSTATIC FORCES IN G-SPACE..."
506 WRITE (iw, '(3f15.9)') ((fgshell_coulomb(j, i), j=1, 3), i=1, SIZE(fgshell_coulomb, 2))
507 WRITE (iw, '(A)') " FIST:: CORE NONBONDED ELECTROSTATIC FORCES IN G-SPACE..."
508 WRITE (iw, '(3f15.9)') ((fgcore_coulomb(j, i), j=1, 3), i=1, SIZE(fgcore_coulomb, 2))
509 END IF
510 END IF
511
512 ! calculate the action of an (external) electric field
513 IF (efield%apply_field) THEN
514 ALLOCATE (ef_f(3, natoms))
515 CALL fist_efield_energy_force(ef_ener, ef_f, ef_pv, atomic_kind_set, particle_set, cell, &
516 efield, use_virial=use_virial, iunit=iw, charges=fist_nonbond_env%charges)
517 END IF
518
519 ! Get intramolecular forces
520 IF (PRESENT(debug)) THEN
521 CALL force_intra_control(molecule_set, molecule_kind_set, &
522 local_molecules, particle_set, shell_particle_set, &
523 core_particle_set, pot_bond, pot_bend, pot_urey_bradley, &
524 pot_torsion, pot_imptors, pot_opbend, pot_shell, pv_bond, pv_bend, &
525 pv_urey_bradley, pv_torsion, pv_imptors, pv_opbend, &
526 debug%f_bond, debug%f_bend, debug%f_torsion, debug%f_ub, &
527 debug%f_imptors, debug%f_opbend, cell, use_virial, atprop_env)
528
529 ELSE
530 CALL force_intra_control(molecule_set, molecule_kind_set, &
531 local_molecules, particle_set, shell_particle_set, &
532 core_particle_set, pot_bond, pot_bend, pot_urey_bradley, &
533 pot_torsion, pot_imptors, pot_opbend, pot_shell, pv_bond, pv_bend, &
534 pv_urey_bradley, pv_torsion, pv_imptors, pv_opbend, &
535 cell=cell, use_virial=use_virial, atprop_env=atprop_env)
536 END IF
537
538 ! Perform global sum of the intra-molecular (bonded) forces for the core-shell atoms
539 IF (shell_present .AND. shell_model_ad) THEN
540 ALLOCATE (fcore_shell_bonded(3, nshell))
541 fcore_shell_bonded = 0.0_dp
542 DO i = 1, natoms
543 shell_index = particle_set(i)%shell_index
544 IF (shell_index /= 0) THEN
545 fcore_shell_bonded(1:3, shell_index) = particle_set(i)%f(1:3)
546 END IF
547 END DO
548 CALL para_env%sum(fcore_shell_bonded)
549 DO i = 1, natoms
550 shell_index = particle_set(i)%shell_index
551 IF (shell_index /= 0) THEN
552 particle_set(i)%f(1:3) = fcore_shell_bonded(1:3, shell_index)
553 END IF
554 END DO
555 DEALLOCATE (fcore_shell_bonded)
556 END IF
557
558 IF (iw > 0) THEN
559 xdum1 = cp_unit_from_cp2k(pot_bond, "kcalmol")
560 xdum2 = cp_unit_from_cp2k(pot_bend, "kcalmol")
561 xdum3 = cp_unit_from_cp2k(pot_urey_bradley, "kcalmol")
562 WRITE (iw, '(A)') " FIST energy contributions in kcal/mol:"
563 WRITE (iw, '(1x,"BOND = ",f13.4,'// &
564 '2x,"ANGLE = ",f13.4,'// &
565 '2x,"UBRAD = ",f13.4)') xdum1, xdum2, xdum3
566 xdum1 = cp_unit_from_cp2k(pot_torsion, "kcalmol")
567 xdum2 = cp_unit_from_cp2k(pot_imptors, "kcalmol")
568 xdum3 = cp_unit_from_cp2k(pot_opbend, "kcalmol")
569 WRITE (iw, '(1x,"TORSION = ",f13.4,'// &
570 '2x,"IMPTORS = ",f13.4,'// &
571 '2x,"OPBEND = ",f13.4)') xdum1, xdum2, xdum3
572
573 WRITE (iw, '(A)') " FIST:: CORRECTED BONDED ELECTROSTATIC FORCES + INTERNAL FORCES..."
574 WRITE (iw, '(3f15.9)') ((particle_set(i)%f(j), j=1, 3), i=1, SIZE(particle_set))
575 IF (shell_present .AND. shell_model_ad) THEN
576 WRITE (iw, '(A)') " FIST:: SHELL CORRECTED BONDED ELECTROSTATIC FORCES + INTERNAL FORCES..."
577 WRITE (iw, '(3f15.9)') ((shell_particle_set(i)%f(j), j=1, 3), i=1, SIZE(shell_particle_set))
578 WRITE (iw, '(A)') " FIST:: CORE CORRECTED BONDED ELECTROSTATIC FORCES + INTERNAL FORCES..."
579 WRITE (iw, '(3f15.9)') ((core_particle_set(i)%f(j), j=1, 3), i=1, SIZE(core_particle_set))
580 END IF
581 END IF
582
583 ! add up all the potential energies
584 thermo%pot = pot_nonbond + pot_bond + pot_bend + pot_torsion + pot_opbend + &
585 pot_imptors + pot_urey_bradley + pot_manybody + pot_shell
586
587 CALL para_env%sum(thermo%pot)
588
589 IF (shell_present) THEN
590 thermo%harm_shell = pot_shell
591 CALL para_env%sum(thermo%harm_shell)
592 END IF
593 ! add g-space contributions if needed
594 IF (ewald_type /= do_ewald_none) THEN
595 ! e_self, e_neut, and ebonded are already summed over all processors
596 ! vg_coulomb is not calculated in parallel
597 thermo%e_gspace = vg_coulomb
598 thermo%pot = thermo%pot + thermo%e_self + thermo%e_neut
599 thermo%pot = thermo%pot + vg_coulomb + thermo%e_bonded
600 ! add the induction energy of the dipoles for polarizable models
601 IF (do_ipol /= do_fist_pol_none) thermo%pot = thermo%pot + thermo%e_induction
602 END IF
603
604 ! add up all the forces
605 !
606 ! nonbonded forces might be calculated for atoms not on this node
607 ! ewald forces are strictly local -> sum only over pnode
608 ! We first sum the forces in f_nonbond, this allows for a more efficient
609 ! global sum in the parallel code and in the end copy them back to part
610 ALLOCATE (f_total(3, natoms))
611 f_total = 0.0_dp
612 DO i = 1, natoms
613 f_total(1, i) = particle_set(i)%f(1) + f_nonbond(1, i)
614 f_total(2, i) = particle_set(i)%f(2) + f_nonbond(2, i)
615 f_total(3, i) = particle_set(i)%f(3) + f_nonbond(3, i)
616 END DO
617 IF (shell_present) THEN
618 ALLOCATE (fshell_total(3, nshell))
619 fshell_total = 0.0_dp
620 ALLOCATE (fcore_total(3, nshell))
621 fcore_total = 0.0_dp
622 DO i = 1, nshell
623 fcore_total(1, i) = core_particle_set(i)%f(1) + fcore_nonbond(1, i)
624 fcore_total(2, i) = core_particle_set(i)%f(2) + fcore_nonbond(2, i)
625 fcore_total(3, i) = core_particle_set(i)%f(3) + fcore_nonbond(3, i)
626 fshell_total(1, i) = shell_particle_set(i)%f(1) + fshell_nonbond(1, i)
627 fshell_total(2, i) = shell_particle_set(i)%f(2) + fshell_nonbond(2, i)
628 fshell_total(3, i) = shell_particle_set(i)%f(3) + fshell_nonbond(3, i)
629 END DO
630 END IF
631
632 IF (iw > 0) THEN
633 WRITE (iw, '(A)') " FIST::(1)INTERNAL + ELECTROSTATIC BONDED + NONBONDED"
634 WRITE (iw, '(3f15.9)') ((f_total(j, i), j=1, 3), i=1, natoms)
635 IF (shell_present .AND. shell_model_ad) THEN
636 WRITE (iw, '(A)') " FIST::(1)SHELL INTERNAL + ELECTROSTATIC BONDED + NONBONDED"
637 WRITE (iw, '(3f15.9)') ((fshell_total(j, i), j=1, 3), i=1, nshell)
638 WRITE (iw, '(A)') " FIST::(1)CORE INTERNAL + ELECTROSTATIC BONDED + NONBONDED"
639 WRITE (iw, '(3f15.9)') ((fcore_total(j, i), j=1, 3), i=1, nshell)
640 END IF
641 END IF
642
643 ! Adding in the reciprocal forces: EWALD is a special case because of distrubted data
644 IF (ewald_type == do_ewald_ewald) THEN
645 node = 0
646 DO iparticle_kind = 1, nparticle_kind
647 nparticle_local = local_particles%n_el(iparticle_kind)
648 DO iparticle_local = 1, nparticle_local
649 ii = local_particles%list(iparticle_kind)%array(iparticle_local)
650 node = node + 1
651 f_total(1, ii) = f_total(1, ii) + fg_coulomb(1, node)
652 f_total(2, ii) = f_total(2, ii) + fg_coulomb(2, node)
653 f_total(3, ii) = f_total(3, ii) + fg_coulomb(3, node)
654 IF (PRESENT(debug)) THEN
655 debug%f_g(1, ii) = debug%f_g(1, ii) + fg_coulomb(1, node)
656 debug%f_g(2, ii) = debug%f_g(2, ii) + fg_coulomb(2, node)
657 debug%f_g(3, ii) = debug%f_g(3, ii) + fg_coulomb(3, node)
658 END IF
659 IF (atprop_env%energy) THEN
660 atprop_env%atener(ii) = atprop_env%atener(ii) + e_coulomb(node)
661 END IF
662 IF (atprop_env%stress) THEN
663 atprop_env%atstress(:, :, ii) = atprop_env%atstress(:, :, ii) + pv_coulomb(:, :, node)
664 END IF
665 END DO
666 END DO
667 IF (atprop_env%energy) THEN
668 DEALLOCATE (e_coulomb)
669 END IF
670 IF (atprop_env%stress) THEN
671 DEALLOCATE (pv_coulomb)
672 END IF
673 END IF
674
675 IF (iw > 0) THEN
676 WRITE (iw, '(A)') " FIST::(2)TOTAL FORCES(1)+ ELECTROSTATIC FORCES"
677 WRITE (iw, '(3f15.9)') ((f_total(j, i), j=1, 3), i=1, natoms)
678 IF (shell_present .AND. shell_model_ad) THEN
679 WRITE (iw, '(A)') " FIST::(2)SHELL TOTAL FORCES(1)+ ELECTROSTATIC FORCES "
680 WRITE (iw, '(3f15.9)') ((fshell_total(j, i), j=1, 3), i=1, nshell)
681 WRITE (iw, '(A)') " FIST::(2)CORE TOTAL FORCES(1)+ ELECTROSTATIC FORCES"
682 WRITE (iw, '(3f15.9)') ((fcore_total(j, i), j=1, 3), i=1, nshell)
683 END IF
684 END IF
685
686 IF (use_virial) THEN
687 ! Add up all the pressure tensors
688 IF (ewald_type == do_ewald_none) THEN
689 virial%pv_virial = pv_nonbond + pv_bond + pv_bend + pv_torsion + pv_imptors + pv_urey_bradley
690 CALL para_env%sum(virial%pv_virial)
691 ELSE
692 ident = 0.0_dp
693 DO i = 1, 3
694 ident(i, i) = 1.0_dp
695 END DO
696
697 virial%pv_virial = pv_nonbond + pv_bond + pv_bend + pv_torsion + pv_imptors + pv_urey_bradley + pv_bc
698 CALL para_env%sum(virial%pv_virial)
699
700 virial%pv_virial = virial%pv_virial + ident*thermo%e_neut
701 virial%pv_virial = virial%pv_virial + pv_g
702 END IF
703 END IF
704
705 ! Sum total forces
706 CALL para_env%sum(f_total)
707 IF (shell_present .AND. shell_model_ad) THEN
708 CALL para_env%sum(fcore_total)
709 CALL para_env%sum(fshell_total)
710 END IF
711
712 ! contributions from fields (currently all quantities are fully replicated!)
713 IF (efield%apply_field) THEN
714 thermo%pot = thermo%pot + ef_ener
715 f_total(1:3, 1:natoms) = f_total(1:3, 1:natoms) + ef_f(1:3, 1:natoms)
716 IF (use_virial) THEN
717 virial%pv_virial(1:3, 1:3) = virial%pv_virial(1:3, 1:3) + ef_pv(1:3, 1:3)
718 END IF
719 DEALLOCATE (ef_f)
720 END IF
721
722 ! Assign to particle_set
723 SELECT CASE (ewald_type)
725 IF (shell_present .AND. shell_model_ad) THEN
726 DO i = 1, natoms
727 shell_index = particle_set(i)%shell_index
728 IF (shell_index == 0) THEN
729 particle_set(i)%f(1:3) = f_total(1:3, i) + fg_coulomb(1:3, i)
730 ELSE
731 atomic_kind => particle_set(i)%atomic_kind
732 CALL get_atomic_kind(atomic_kind=atomic_kind, shell=shell, mass=mass)
733 fc = shell%mass_core/mass
734 fcore_total(1:3, shell_index) = fcore_total(1:3, shell_index) + particle_set(i)%f(1:3)*fc
735 fs = shell%mass_shell/mass
736 fshell_total(1:3, shell_index) = fshell_total(1:3, shell_index) + particle_set(i)%f(1:3)*fs
737 END IF
738 END DO
739
740 DO i = 1, nshell
741 core_particle_set(i)%f(1) = fcore_total(1, i) + fgcore_coulomb(1, i)
742 core_particle_set(i)%f(2) = fcore_total(2, i) + fgcore_coulomb(2, i)
743 core_particle_set(i)%f(3) = fcore_total(3, i) + fgcore_coulomb(3, i)
744 shell_particle_set(i)%f(1) = fshell_total(1, i) + fgshell_coulomb(1, i)
745 shell_particle_set(i)%f(2) = fshell_total(2, i) + fgshell_coulomb(2, i)
746 shell_particle_set(i)%f(3) = fshell_total(3, i) + fgshell_coulomb(3, i)
747 END DO
748
749 ELSEIF (shell_present .AND. .NOT. shell_model_ad) THEN
750 cpabort("Non adiabatic shell-model not implemented.")
751 ELSE
752 DO i = 1, natoms
753 particle_set(i)%f(1) = f_total(1, i) + fg_coulomb(1, i)
754 particle_set(i)%f(2) = f_total(2, i) + fg_coulomb(2, i)
755 particle_set(i)%f(3) = f_total(3, i) + fg_coulomb(3, i)
756 END DO
757 END IF
759 IF (shell_present .AND. shell_model_ad) THEN
760 DO i = 1, natoms
761 shell_index = particle_set(i)%shell_index
762 IF (shell_index == 0) THEN
763 particle_set(i)%f(1:3) = f_total(1:3, i)
764 ELSE
765 atomic_kind => particle_set(i)%atomic_kind
766 CALL get_atomic_kind(atomic_kind=atomic_kind, shell=shell, mass=mass)
767 fc = shell%mass_core/mass
768 fcore_total(1:3, shell_index) = fcore_total(1:3, shell_index) + particle_set(i)%f(1:3)*fc
769 fs = shell%mass_shell/mass
770 fshell_total(1:3, shell_index) = fshell_total(1:3, shell_index) + particle_set(i)%f(1:3)*fs
771 END IF
772 END DO
773 DO i = 1, nshell
774 core_particle_set(i)%f(1) = fcore_total(1, i)
775 core_particle_set(i)%f(2) = fcore_total(2, i)
776 core_particle_set(i)%f(3) = fcore_total(3, i)
777 shell_particle_set(i)%f(1) = fshell_total(1, i)
778 shell_particle_set(i)%f(2) = fshell_total(2, i)
779 shell_particle_set(i)%f(3) = fshell_total(3, i)
780 END DO
781 ELSEIF (shell_present .AND. .NOT. shell_model_ad) THEN
782 cpabort("Non adiabatic shell-model not implemented.")
783 ELSE
784 DO i = 1, natoms
785 particle_set(i)%f(1) = f_total(1, i)
786 particle_set(i)%f(2) = f_total(2, i)
787 particle_set(i)%f(3) = f_total(3, i)
788 END DO
789 END IF
790 END SELECT
791
792 IF (iw > 0) THEN
793 WRITE (iw, '(A)') " FIST::(3)TOTAL FORCES - THE END..."
794 WRITE (iw, '(3f15.9)') ((particle_set(i)%f(j), j=1, 3), i=1, natoms)
795 IF (shell_present .AND. shell_model_ad) THEN
796 WRITE (iw, '(A)') " FIST::(3)SHELL TOTAL FORCES - THE END..."
797 WRITE (iw, '(3f15.9)') ((shell_particle_set(i)%f(j), j=1, 3), i=1, nshell)
798 WRITE (iw, '(A)') " FIST::(3)CORE TOTAL FORCES - THE END..."
799 WRITE (iw, '(3f15.9)') ((core_particle_set(i)%f(j), j=1, 3), i=1, nshell)
800 END IF
801 WRITE (iw, '(A,f15.9)') "Energy after FIST calculation.. exiting now ::", thermo%pot
802 END IF
803 !
804 ! If we are doing debugging, check if variables are present and assign
805 !
806 IF (PRESENT(debug)) THEN
807 CALL para_env%sum(pot_manybody)
808 debug%pot_manybody = pot_manybody
809 CALL para_env%sum(pot_nonbond)
810 debug%pot_nonbond = pot_nonbond
811 CALL para_env%sum(pot_bond)
812 debug%pot_bond = pot_bond
813 CALL para_env%sum(pot_bend)
814 debug%pot_bend = pot_bend
815 CALL para_env%sum(pot_torsion)
816 debug%pot_torsion = pot_torsion
817 CALL para_env%sum(pot_imptors)
818 debug%pot_imptors = pot_imptors
819 CALL para_env%sum(pot_opbend)
820 debug%pot_opbend = pot_opbend
821 CALL para_env%sum(pot_urey_bradley)
822 debug%pot_urey_bradley = pot_urey_bradley
823 CALL para_env%sum(f_nonbond)
824 debug%f_nonbond = f_nonbond
825 IF (use_virial) THEN
826 CALL para_env%sum(pv_nonbond)
827 debug%pv_nonbond = pv_nonbond
828 CALL para_env%sum(pv_bond)
829 debug%pv_bond = pv_bond
830 CALL para_env%sum(pv_bend)
831 debug%pv_bend = pv_bend
832 CALL para_env%sum(pv_torsion)
833 debug%pv_torsion = pv_torsion
834 CALL para_env%sum(pv_imptors)
835 debug%pv_imptors = pv_imptors
836 CALL para_env%sum(pv_urey_bradley)
837 debug%pv_ub = pv_urey_bradley
838 END IF
839 SELECT CASE (ewald_type)
841 debug%pot_g = vg_coulomb
842 debug%pv_g = pv_g
843 debug%f_g = fg_coulomb
844 CASE (do_ewald_ewald)
845 debug%pot_g = vg_coulomb
846 IF (use_virial) debug%pv_g = pv_g
847 CASE default
848 debug%pot_g = 0.0_dp
849 debug%f_g = 0.0_dp
850 IF (use_virial) debug%pv_g = 0.0_dp
851 END SELECT
852 END IF
853
854 ! print properties if requested
855 print_section => section_vals_get_subs_vals(mm_section, "PRINT")
856 CALL print_fist(fist_env, print_section, atomic_kind_set, particle_set, cell)
857
858 ! deallocating all local variables
859 IF (ALLOCATED(fg_coulomb)) THEN
860 DEALLOCATE (fg_coulomb)
861 END IF
862 IF (ALLOCATED(f_total)) THEN
863 DEALLOCATE (f_total)
864 END IF
865 DEALLOCATE (f_nonbond)
866 IF (shell_present) THEN
867 DEALLOCATE (fshell_total)
868 IF (ewald_type /= do_ewald_none) THEN
869 DEALLOCATE (fgshell_coulomb)
870 END IF
871 DEALLOCATE (fshell_nonbond)
872 END IF
873 CALL cp_print_key_finished_output(iw, logger, mm_section, "PRINT%DERIVATIVES")
874 CALL cp_print_key_finished_output(iw2, logger, mm_section, "PRINT%EWALD_INFO")
875 CALL timestop(handle)
876
877 END SUBROUTINE fist_calc_energy_force
878
879! **************************************************************************************************
880!> \brief Print properties number according the requests in input file
881!> \param fist_env ...
882!> \param print_section ...
883!> \param atomic_kind_set ...
884!> \param particle_set ...
885!> \param cell ...
886!> \par History
887!> [01.2006] created
888!> \author Teodoro Laino
889! **************************************************************************************************
890 SUBROUTINE print_fist(fist_env, print_section, atomic_kind_set, particle_set, cell)
891 TYPE(fist_environment_type), POINTER :: fist_env
892 TYPE(section_vals_type), POINTER :: print_section
893 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
894 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
895 TYPE(cell_type), POINTER :: cell
896
897 INTEGER :: unit_nr
898 TYPE(cp_logger_type), POINTER :: logger
899 TYPE(fist_nonbond_env_type), POINTER :: fist_nonbond_env
900 TYPE(section_vals_type), POINTER :: print_key
901
902 NULLIFY (logger, print_key, fist_nonbond_env)
903 logger => cp_get_default_logger()
904 print_key => section_vals_get_subs_vals(print_section, "dipole")
905 CALL fist_env_get(fist_env, fist_nonbond_env=fist_nonbond_env)
906 IF (btest(cp_print_key_should_output(logger%iter_info, print_key), &
907 cp_p_file)) THEN
908 unit_nr = cp_print_key_unit_nr(logger, print_section, "dipole", &
909 extension=".data", middle_name="MM_DIPOLE", log_filename=.false.)
910 CALL fist_dipole(fist_env, print_section, atomic_kind_set, particle_set, &
911 cell, unit_nr, fist_nonbond_env%charges)
912 CALL cp_print_key_finished_output(unit_nr, logger, print_key)
913 END IF
914
915 END SUBROUTINE print_fist
916
917END MODULE fist_force
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.
Holds information on atomic properties.
Handles all functions related to the CELL.
subroutine, public init_cell(cell, hmat, periodic)
Initialise/readjust a simulation cell after hmat has been changed.
Handles all functions related to the CELL.
Definition cell_types.F:15
various routines to log and control the output. The idea is that decisions about where to log should ...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
integer, parameter, public cp_p_file
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
types that represent a subsys, i.e. a part of the system
subroutine, public cp_subsys_get(subsys, ref_count, atomic_kinds, atomic_kind_set, particles, particle_set, local_particles, molecules, molecule_set, molecule_kinds, molecule_kind_set, local_molecules, para_env, colvar_p, shell_particles, core_particles, gci, multipoles, natom, nparticle, ncore, nshell, nkind, atprop, virial, results, cell)
returns information about various attributes of the given subsys
unit conversion facility
Definition cp_units.F:30
real(kind=dp) function, public cp_unit_from_cp2k(value, unit_str, defaults, power)
converts from the internal cp2k units to the given unit
Definition cp_units.F:1179
stores a lists of integer that are local to a processor. The idea is that these integers represent ob...
subroutine, public ewald_env_get(ewald_env, ewald_type, alpha, eps_pol, epsilon, gmax, ns_max, o_spline, group, para_env, poisson_section, precs, rcut, do_multipoles, max_multipole, do_ipol, max_ipol_iter, interaction_cutoffs, cell_hmat)
Purpose: Get the EWALD environment.
subroutine, public ewald_pw_grid_update(ewald_pw, ewald_env, cell_hmat)
Rescales pw_grids for given box, if necessary.
Treats the electrostatic for multipoles (up to quadrupoles)
recursive subroutine, public ewald_multipole_evaluate(ewald_env, ewald_pw, nonbond_env, cell, particle_set, local_particles, energy_local, energy_glob, e_neut, e_self, task, do_correction_bonded, do_forces, do_stress, do_efield, radii, charges, dipoles, quadrupoles, forces_local, forces_glob, pv_local, pv_glob, efield0, efield1, efield2, iw, do_debug, atomic_kind_set, mm_section)
Computes the potential and the force for a lattice sum of multipoles (up to quadrupole)
subroutine, public ewald_evaluate(ewald_env, ewald_pw, cell, atomic_kind_set, particle_set, local_particles, fg_coulomb, vg_coulomb, pv_g, use_virial, charges, e_coulomb, pv_coulomb)
computes the potential and the force from the g-space part of the 1/r potential Ref....
Definition ewalds.F:79
subroutine, public ewald_self_atom(ewald_env, atomic_kind_set, local_particles, e_self, charges)
Computes the self interaction per atom.
Definition ewalds.F:390
subroutine, public ewald_print(iw, pot_nonbond, e_gspace, e_self, e_neut, e_bonded)
...
Definition ewalds.F:450
subroutine, public ewald_self(ewald_env, cell, atomic_kind_set, local_particles, e_self, e_neut, charges)
Computes the self interaction from g-space and the neutralizing background.
Definition ewalds.F:305
an exclusion type
subroutine, public fist_efield_energy_force(qenergy, qforce, qpv, atomic_kind_set, particle_set, cell, efield, use_virial, iunit, charges)
...
subroutine, public fist_dipole(fist_env, print_section, atomic_kind_set, particle_set, cell, unit_nr, charges)
Evaluates the Dipole of a classical charge distribution(point-like) possibly using the berry phase fo...
subroutine, public fist_env_get(fist_env, atomic_kind_set, particle_set, ewald_pw, local_particles, local_molecules, molecule_kind_set, molecule_set, cell, cell_ref, ewald_env, fist_nonbond_env, thermo, para_env, subsys, qmmm, qmmm_env, input, shell_model, shell_model_ad, shell_particle_set, core_particle_set, multipoles, results, exclusions, efield)
Purpose: Get the FIST environment.
subroutine, public fist_calc_energy_force(fist_env, debug)
Calculates the total potential energy, total force, and the total pressure tensor from the potentials...
Definition fist_force.F:112
subroutine, public force_intra_control(molecule_set, molecule_kind_set, local_molecules, particle_set, shell_particle_set, core_particle_set, pot_bond, pot_bend, pot_urey_bradley, pot_torsion, pot_imp_torsion, pot_opbend, pot_shell, pv_bond, pv_bend, pv_urey_bradley, pv_torsion, pv_imp_torsion, pv_opbend, f_bond, f_bend, f_torsion, f_ub, f_imptor, f_opbend, cell, use_virial, atprop_env)
Computes the the intramolecular energies, forces, and pressure tensors.
subroutine, public list_control(atomic_kind_set, particle_set, local_particles, cell, fist_nonbond_env, para_env, mm_section, shell_particle_set, core_particle_set, force_update, exclusions)
...
subroutine, public force_nonbond(fist_nonbond_env, ewald_env, particle_set, cell, pot_nonbond, f_nonbond, pv_nonbond, fshell_nonbond, fcore_nonbond, atprop_env, atomic_kind_set, use_virial)
Calculates the force and the potential of the minimum image, and the pressure tensor.
subroutine, public bonded_correct_gaussian(fist_nonbond_env, atomic_kind_set, local_particles, particle_set, ewald_env, v_bonded_corr, pv_bc, shell_particle_set, core_particle_set, atprop_env, cell, use_virial)
corrects electrostatics for bonded terms
subroutine, public fist_pol_evaluate(atomic_kind_set, multipoles, ewald_env, ewald_pw, fist_nonbond_env, cell, particle_set, local_particles, thermo, vg_coulomb, pot_nonbond, f_nonbond, fg_coulomb, use_virial, pv_g, pv_nonbond, mm_section, do_ipol)
Main driver for evaluating energy/forces in a polarizable forcefield.
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public do_fist_pol_none
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
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
subroutine, public density_nonbond(fist_nonbond_env, particle_set, cell, para_env)
...
subroutine, public energy_manybody(fist_nonbond_env, atomic_kind_set, local_particles, particle_set, cell, pot_manybody, para_env, mm_section)
computes the embedding contribution to the energy
subroutine, public force_nonbond_manybody(fist_nonbond_env, particle_set, cell, f_nonbond, pv_nonbond, use_virial)
...
Interface to the message passing library MPI.
Define the molecule kind structure types and the corresponding functionality.
Define the data structure for the molecule information.
Multipole structure: for multipole (fixed and induced) in FF based MD.
Define the data structure for the particle information.
Definition pme.F:13
subroutine, public pme_evaluate(ewald_env, ewald_pw, box, particle_set, vg_coulomb, fg_coulomb, pv_g, shell_particle_set, core_particle_set, fgshell_coulomb, fgcore_coulomb, use_virial, charges, atprop)
...
Definition pme.F:98
functions related to the poisson solver on regular grids
integer, parameter, public do_ewald_pme
integer, parameter, public do_ewald_ewald
integer, parameter, public do_ewald_none
integer, parameter, public do_ewald_spme
Calculate the electrostatic energy by the Smooth Particle Ewald method.
Definition spme.F:14
subroutine, public spme_evaluate(ewald_env, ewald_pw, box, particle_set, fg_coulomb, vg_coulomb, pv_g, shell_particle_set, core_particle_set, fgshell_coulomb, fgcore_coulomb, use_virial, charges, atprop)
...
Definition spme.F:95
Provides all information about an atomic kind.
type for the atomic properties
Type defining parameters related to the simulation cell.
Definition cell_types.F:55
type of a logger, at the moment it contains just a print level starting at which level it should be l...
represents a system: atoms, molecules, their pos,vel,...
structure to store local (to a processor) ordered lists of integers.
A type used to store lists of exclusions and onfos.
stores all the informations relevant to an mpi environment