53 #include "./base/base_uses.f90"
58 LOGICAL,
PRIVATE :: debug_this_module = .false.
59 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qmmmx_util'
75 TYPE(qmmmx_env_type),
POINTER :: qmmmx_env
78 TYPE(cell_type),
POINTER :: cell_core, cell_extended
79 TYPE(cp_subsys_type),
POINTER :: subsys_core, subsys_extended
80 TYPE(particle_type),
DIMENSION(:),
POINTER :: particles_core, particles_extended
82 NULLIFY (cell_core, cell_extended)
83 NULLIFY (subsys_core, subsys_extended)
84 NULLIFY (particles_core, particles_extended)
94 particles_extended => subsys_extended%particles%els
95 particles_core => subsys_core%particles%els
96 DO ip = 1,
SIZE(particles_extended)
97 particles_core(ip)%r = particles_extended(ip)%r
116 TYPE(cp_subsys_type),
POINTER :: subsys
117 TYPE(section_vals_type),
POINTER :: qmmm_section
118 LOGICAL,
OPTIONAL :: labels_changed
120 CHARACTER(LEN=default_string_length),
POINTER :: adaptive_exclude_molecules(:)
121 INTEGER :: i_rep_section, i_rep_val, ip, max_n_qm, n_new, n_rep_exclude, n_rep_section, &
122 n_rep_val, natoms, output_unit, qm_extended_seed_min_label_val
123 INTEGER,
ALLOCATABLE :: new_full_labels(:), orig_full_labels(:)
124 INTEGER,
POINTER :: broken_bonds(:), cur_indices(:), &
125 cur_labels(:), mm_index_entry(:), &
126 new_indices(:), new_labels(:)
127 LOGICAL :: explicit, qm_extended_seed_is_core_list
128 REAL(
dp),
ALLOCATABLE :: nearest_dist(:)
129 REAL(
dp),
POINTER :: r_buf(:), r_core(:), r_qm(:)
130 TYPE(cell_type),
POINTER :: cell
131 TYPE(fist_neighbor_type),
POINTER :: nlist
132 TYPE(molecule_list_type),
POINTER :: molecules
133 TYPE(molecule_type),
DIMENSION(:),
POINTER :: molecule_set
134 TYPE(mp_para_env_type),
POINTER :: para_env
135 TYPE(particle_list_type),
POINTER :: particles
136 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
137 TYPE(section_vals_type),
POINTER :: force_mixing_section, &
138 non_adaptive_section, qm_kind_section, &
143 IF (debug_this_module .AND. output_unit > 0)
WRITE (output_unit, *)
"BOB starting update_force_mixing_labels"
146 CALL get_force_mixing_indices(force_mixing_section, cur_indices, cur_labels)
147 IF (debug_this_module .AND. output_unit > 0)
WRITE (output_unit, *)
"BOB got cur_indices ",
SIZE(cur_indices)
148 IF (debug_this_module .AND. output_unit > 0)
WRITE (output_unit, *)
"BOB got cur_labels ",
SIZE(cur_labels)
152 NULLIFY (r_core, r_qm, r_buf, adaptive_exclude_molecules, broken_bonds)
156 l_val=qm_extended_seed_is_core_list)
160 CALL section_vals_val_get(force_mixing_section,
"ADAPTIVE_EXCLUDE_MOLECULES", n_rep_val=n_rep_exclude)
161 IF (n_rep_exclude > 0)
THEN
162 CALL section_vals_val_get(force_mixing_section,
"ADAPTIVE_EXCLUDE_MOLECULES", c_vals=adaptive_exclude_molecules)
170 NULLIFY (particles, molecules)
171 CALL cp_subsys_get(subsys=subsys, particles=particles, molecules=molecules)
172 particle_set => particles%els
173 molecule_set => molecules%els
175 natoms =
SIZE(particle_set)
178 NULLIFY (new_indices, new_labels)
179 CALL reallocate(new_indices, 1,
SIZE(cur_indices))
180 CALL reallocate(new_labels, 1,
SIZE(cur_labels))
183 ALLOCATE (new_full_labels(natoms))
190 CALL make_neighbor_list(force_mixing_section, subsys, cell, max(r_core(2), r_qm(2), r_buf(2)), nlist)
193 NULLIFY (mm_index_entry)
197 DO i_rep_section = 1, n_rep_section
198 CALL section_vals_val_get(qm_kind_section,
"MM_INDEX", i_rep_section=i_rep_section, n_rep_val=n_rep_val)
199 DO i_rep_val = 1, n_rep_val
200 CALL section_vals_val_get(qm_kind_section,
"MM_INDEX", i_rep_section=i_rep_section, i_rep_val=i_rep_val, &
201 i_vals=mm_index_entry)
202 DO ip = 1,
SIZE(mm_index_entry)
204 new_full_labels, max_n_qm)
209 IF (debug_this_module .AND. output_unit > 0)
THEN
210 WRITE (output_unit, *)
"BOB core_list new_indices ", new_indices(1:n_new)
211 WRITE (output_unit, *)
"BOB core_list new_labels ", new_labels(1:n_new)
216 can_return_null=.true.)
217 IF (
ASSOCIATED(non_adaptive_section))
THEN
220 DO i_rep_section = 1, n_rep_section
221 CALL section_vals_val_get(qm_kind_section,
"MM_INDEX", i_rep_section=i_rep_section, n_rep_val=n_rep_val)
222 DO i_rep_val = 1, n_rep_val
223 CALL section_vals_val_get(qm_kind_section,
"MM_INDEX", i_rep_section=i_rep_section, i_rep_val=i_rep_val, &
224 i_vals=mm_index_entry)
225 DO ip = 1,
SIZE(mm_index_entry)
227 new_full_labels, max_n_qm)
232 IF (debug_this_module .AND. output_unit > 0)
THEN
233 WRITE (output_unit, *)
"BOB core_list + non adaptive QM new_indices ", new_indices(1:n_new)
234 WRITE (output_unit, *)
"BOB core_list + non adaptive QM new_labels ", new_labels(1:n_new)
237 can_return_null=.true.)
238 IF (
ASSOCIATED(non_adaptive_section))
THEN
241 DO i_rep_section = 1, n_rep_section
242 CALL section_vals_val_get(qm_kind_section,
"MM_INDEX", i_rep_section=i_rep_section, n_rep_val=n_rep_val)
243 DO i_rep_val = 1, n_rep_val
244 CALL section_vals_val_get(qm_kind_section,
"MM_INDEX", i_rep_section=i_rep_section, i_rep_val=i_rep_val, &
245 i_vals=mm_index_entry)
246 DO ip = 1,
SIZE(mm_index_entry)
248 new_full_labels, max_n_qm)
254 IF (debug_this_module .AND. output_unit > 0)
THEN
255 WRITE (output_unit, *)
"BOB core_list + non adaptive QM+buffer new_indices ", new_indices(1:n_new)
256 WRITE (output_unit, *)
"BOB core_list + non adaptive QM+buffer new_labels ", new_labels(1:n_new)
260 ALLOCATE (nearest_dist(natoms))
263 ALLOCATE (orig_full_labels(natoms))
265 orig_full_labels(cur_indices(:)) = cur_labels(:)
272 CALL add_layer_hysteretically( &
273 nlist, particle_set, cell, nearest_dist, &
274 orig_full_labels, new_full_labels, n_new, new_indices, new_labels, &
276 max_n_qm, adaptive_exclude_molecules, molecule_set, broken_bonds)
278 DEALLOCATE (broken_bonds)
280 IF (debug_this_module .AND. output_unit > 0)
THEN
281 WRITE (output_unit, *)
"BOB core new_indices ", new_indices(1:n_new)
282 WRITE (output_unit, *)
"BOB core new_labels ", new_labels(1:n_new)
288 IF (debug_this_module .AND. output_unit > 0) &
289 WRITE (output_unit, *)
"BOB QM_extended_seed_is_core_list ", qm_extended_seed_is_core_list
290 IF (qm_extended_seed_is_core_list)
THEN
295 CALL add_layer_hysteretically(nlist, particle_set, cell, nearest_dist, &
296 orig_full_labels, new_full_labels, n_new, new_indices, new_labels, &
299 max_n_qm, adaptive_exclude_molecules, molecule_set)
301 IF (debug_this_module .AND. output_unit > 0)
THEN
302 WRITE (output_unit, *)
"BOB extended new_indices ", new_indices(1:n_new)
303 WRITE (output_unit, *)
"BOB extended new_labels ", new_labels(1:n_new)
307 CALL add_layer_hysteretically( &
308 nlist, particle_set, cell, nearest_dist, &
309 orig_full_labels, new_full_labels, n_new, new_indices, new_labels, &
311 max_n_qm, adaptive_exclude_molecules, molecule_set, broken_bonds)
313 DEALLOCATE (broken_bonds)
315 IF (debug_this_module .AND. output_unit > 0)
THEN
316 WRITE (output_unit, *)
"BOB buffer new_indices ", new_indices(1:n_new)
317 WRITE (output_unit, *)
"BOB buffer new_labels ", new_labels(1:n_new)
320 DEALLOCATE (nearest_dist)
322 IF (
PRESENT(labels_changed)) labels_changed = any(new_full_labels /= orig_full_labels)
324 DEALLOCATE (orig_full_labels)
325 DEALLOCATE (new_full_labels)
328 CALL reallocate(new_indices, 1, n_new)
329 CALL reallocate(new_labels, 1, n_new)
338 DEALLOCATE (cur_indices, cur_labels)
343 IF (para_env%is_source() .AND. output_unit > 0)
THEN
344 WRITE (unit=output_unit, fmt=
'(A,A,I6,A,I5,A,I5,A,I5)') &
345 "QMMM FORCE MIXING final count (not including links): ", &
366 SUBROUTINE add_new_label(ip, label, n_new, new_indices, new_labels, new_full_labels, max_n_qm)
367 INTEGER :: ip, label, n_new
368 INTEGER,
POINTER :: new_indices(:), new_labels(:)
369 INTEGER :: new_full_labels(:), max_n_qm
371 INTEGER :: i, old_index
376 IF (new_indices(i) == ip)
THEN
381 IF (old_index <= 0) &
382 CALL cp_abort(__location__, &
383 "add_new_label found atom with a label "// &
384 "already set, but not in new_indices array")
385 new_labels(old_index) = label
388 IF (n_new > max_n_qm) &
389 CALL cp_abort(__location__, &
390 "add_new_label tried to add more atoms "// &
391 "than allowed by &FORCE_MIXING&MAX_N_QM!")
392 IF (n_new >
SIZE(new_indices))
CALL reallocate(new_indices, 1, n_new + 9)
393 IF (n_new >
SIZE(new_labels))
CALL reallocate(new_labels, 1, n_new + 9)
394 new_indices(n_new) = ip
395 new_labels(n_new) = label
397 new_full_labels(ip) = label
398 END SUBROUTINE add_new_label
420 SUBROUTINE add_layer_hysteretically(nlist, particle_set, cell, nearest_dist, &
421 orig_full_labels, new_full_labels, n_new, new_indices, new_labels, &
422 seed_min_label_val, seed_max_label_val, set_label_val, r_inout, max_n_qm, &
423 adaptive_exclude_molecules, molecule_set, broken_bonds)
424 TYPE(fist_neighbor_type),
POINTER :: nlist
425 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
426 TYPE(cell_type),
POINTER :: cell
427 REAL(
dp) :: nearest_dist(:)
428 INTEGER :: orig_full_labels(:), new_full_labels(:), &
430 INTEGER,
POINTER :: new_indices(:), new_labels(:)
431 INTEGER :: seed_min_label_val, seed_max_label_val, &
433 REAL(
dp) :: r_inout(2)
435 CHARACTER(len=*),
POINTER :: adaptive_exclude_molecules(:)
436 TYPE(molecule_type),
DIMENSION(:),
OPTIONAL, &
437 POINTER :: molecule_set
438 INTEGER,
OPTIONAL,
POINTER :: broken_bonds(:)
440 INTEGER :: i_ind, im, im_exclude, ip, ipair, &
441 ipairkind, j_ind, output_unit
442 LOGICAL :: adaptive_exclude, i_in_new_seed, i_outside_new_seed, j_in_new_seed, &
443 j_outside_new_seed, molec_in_inner, molec_in_outer
444 REAL(
dp) :: r_ij(3), r_ij_mag
448 IF (debug_this_module .AND. output_unit > 0)
WRITE (output_unit, *)
"BOB adding hysteretically seed ", &
449 seed_min_label_val, seed_max_label_val,
" set ", set_label_val,
" r ", r_inout
451 nearest_dist = huge(1.0_dp)
453 DO ipairkind = 1,
SIZE(nlist%neighbor_kind_pairs)
454 DO ipair = 1, nlist%neighbor_kind_pairs(ipairkind)%npairs
456 i_ind = nlist%neighbor_kind_pairs(ipairkind)%list(1, ipair)
457 j_ind = nlist%neighbor_kind_pairs(ipairkind)%list(2, ipair)
459 i_in_new_seed = (new_full_labels(i_ind) >= seed_min_label_val .AND. new_full_labels(i_ind) <= seed_max_label_val)
460 i_outside_new_seed = (new_full_labels(i_ind) < seed_min_label_val)
461 j_in_new_seed = (new_full_labels(j_ind) >= seed_min_label_val .AND. new_full_labels(j_ind) <= seed_max_label_val)
462 j_outside_new_seed = (new_full_labels(j_ind) < seed_min_label_val)
464 IF ((i_in_new_seed .AND. j_outside_new_seed) .OR. (j_in_new_seed .AND. i_outside_new_seed))
THEN
465 r_ij =
pbc(particle_set(i_ind)%r - particle_set(j_ind)%r, cell)
466 r_ij_mag = sqrt(sum(r_ij**2))
467 IF (i_in_new_seed .AND. j_outside_new_seed .AND. (r_ij_mag < nearest_dist(j_ind)))
THEN
468 nearest_dist(j_ind) = r_ij_mag
470 IF (j_in_new_seed .AND. i_outside_new_seed .AND. (r_ij_mag < nearest_dist(i_ind)))
THEN
471 nearest_dist(i_ind) = r_ij_mag
480 DO im = 1,
SIZE(molecule_set)
482 IF (
ASSOCIATED(adaptive_exclude_molecules))
THEN
483 adaptive_exclude = .false.
484 DO im_exclude = 1,
SIZE(adaptive_exclude_molecules)
485 IF (trim(molecule_set(im)%molecule_kind%name) == trim(adaptive_exclude_molecules(im_exclude)) .OR. &
486 trim(molecule_set(im)%molecule_kind%name) ==
'_QM_'//trim(adaptive_exclude_molecules(im_exclude))) &
487 adaptive_exclude = .true.
489 IF (adaptive_exclude) cycle
491 molec_in_inner = any(nearest_dist(molecule_set(im)%first_atom:molecule_set(im)%last_atom) <= r_inout(1))
492 molec_in_outer = any(nearest_dist(molecule_set(im)%first_atom:molecule_set(im)%last_atom) <= r_inout(2))
493 IF (molec_in_inner)
THEN
494 DO ip = molecule_set(im)%first_atom, molecule_set(im)%last_atom
496 IF (new_full_labels(ip) < set_label_val) &
497 CALL add_new_label(ip, set_label_val, n_new, new_indices, new_labels, new_full_labels, max_n_qm)
499 ELSE IF (molec_in_outer)
THEN
500 IF (any(orig_full_labels(molecule_set(im)%first_atom:molecule_set(im)%last_atom) >= set_label_val))
THEN
501 DO ip = molecule_set(im)%first_atom, molecule_set(im)%last_atom
503 IF (new_full_labels(ip) < set_label_val) &
504 CALL add_new_label(ip, set_label_val, n_new, new_indices, new_labels, new_full_labels, max_n_qm)
509 IF (
PRESENT(broken_bonds))
CALL reallocate(broken_bonds, 1, 0)
511 END SUBROUTINE add_layer_hysteretically
521 SUBROUTINE make_neighbor_list(force_mixing_section, subsys, cell, r_max, nlist)
522 TYPE(section_vals_type),
POINTER :: force_mixing_section
523 TYPE(cp_subsys_type),
POINTER :: subsys
524 TYPE(cell_type),
POINTER :: cell
526 TYPE(fist_neighbor_type),
POINTER :: nlist
528 CHARACTER(LEN=default_string_length) :: kind_name
529 CHARACTER(LEN=default_string_length),
POINTER :: kind_name_a(:)
532 REAL(
dp),
ALLOCATABLE :: r_max_a(:, :), r_minsq_a(:, :)
533 TYPE(atomic_kind_type),
POINTER :: atomic_kind
535 ALLOCATE (r_max_a(
SIZE(subsys%atomic_kinds%els),
SIZE(subsys%atomic_kinds%els)))
536 ALLOCATE (r_minsq_a(
SIZE(subsys%atomic_kinds%els),
SIZE(subsys%atomic_kinds%els)))
538 r_minsq_a = epsilon(1.0_dp)
541 ALLOCATE (kind_name_a(
SIZE(subsys%atomic_kinds%els)))
542 DO ik = 1,
SIZE(subsys%atomic_kinds%els)
543 atomic_kind => subsys%atomic_kinds%els(ik)
545 kind_name_a(ik) = kind_name
550 DO ik = 1,
SIZE(subsys%atomic_kinds%els)
551 atomic_kind => subsys%atomic_kinds%els(ik)
561 cell=cell, r_max=r_max_a, r_minsq=r_minsq_a, &
562 ei_scale14=1.0_dp, vdw_scale14=1.0_dp, nonbonded=nlist, &
563 para_env=subsys%para_env, build_from_scratch=.true., geo_check=.false., &
564 mm_section=force_mixing_section)
566 DEALLOCATE (r_max_a, r_minsq_a)
569 DO ik = 1,
SIZE(subsys%atomic_kinds%els)
572 DEALLOCATE (kind_name_a)
574 END SUBROUTINE make_neighbor_list
587 TYPE(cp_subsys_type),
POINTER :: subsys
588 TYPE(section_vals_type),
POINTER :: qmmm_section, qmmm_core_section, &
589 qmmm_extended_section
591 CHARACTER(len=default_string_length),
POINTER :: elem_mapping(:, :), elem_mapping_entry(:)
592 INTEGER :: delta_charge, i_rep_section_core, i_rep_section_extended, i_rep_val_core, &
593 i_rep_val_extended, ielem, ip, n_elements, output_unit
594 INTEGER,
POINTER :: cur_indices(:), cur_labels(:)
595 LOGICAL :: mapped, new_element_core, &
597 TYPE(particle_type),
DIMENSION(:),
POINTER :: particles
598 TYPE(section_vals_type),
POINTER :: buffer_non_adaptive_section, &
600 force_mixing_section, link_section, &
603 NULLIFY (qmmm_core_section, qmmm_extended_section)
612 IF (
ASSOCIATED(link_section))
THEN
617 can_return_null=.true.)
619 IF (
ASSOCIATED(link_section))
THEN
620 NULLIFY (dup_link_section)
626 IF (debug_this_module .AND. output_unit > 0)
THEN
628 WRITE (output_unit, *)
"core section has LINKs ",
ASSOCIATED(link_section)
631 WRITE (output_unit, *)
"extended section has LINKs ",
ASSOCIATED(link_section)
638 CALL section_vals_val_get(force_mixing_section,
"QM_KIND_ELEMENT_MAPPING", n_rep_val=n_elements)
639 ALLOCATE (elem_mapping(2, n_elements))
640 DO ielem = 1, n_elements
641 CALL section_vals_val_get(force_mixing_section,
"QM_KIND_ELEMENT_MAPPING", i_rep_val=ielem, c_vals=elem_mapping_entry)
642 elem_mapping(1:2, ielem) = elem_mapping_entry(1:2)
646 CALL get_force_mixing_indices(force_mixing_section, cur_indices, cur_labels)
647 IF (
SIZE(cur_indices) <= 0) &
648 cpabort(
"cur_indices is empty, found no QM atoms")
650 IF (debug_this_module .AND. output_unit > 0)
THEN
651 WRITE (output_unit, *)
"cur_indices ", cur_indices
652 WRITE (output_unit, *)
"cur_labels ", cur_labels
656 particles => subsys%particles%els
658 DO ip = 1,
SIZE(cur_indices)
662 DO ielem = 1, n_elements
663 IF (trim(particles(cur_indices(ip))%atomic_kind%element_symbol) == trim(elem_mapping(1, ielem)))
THEN
669 CALL cp_abort(__location__, &
670 "Force-mixing failed to find QM_KIND mapping for atom of type "// &
671 trim(particles(cur_indices(ip))%atomic_kind%element_symbol)// &
679 IF (i_rep_section_core <= 0) &
680 CALL cp_abort(__location__, &
681 "Force-mixing QM didn't find any QM_KIND sections, "// &
682 "so no core specified!")
683 i_rep_section_extended = i_rep_section_core
684 DO ielem = 1, n_elements
685 new_element_core = .true.
686 new_element_extended = .true.
687 DO ip = 1,
SIZE(cur_indices)
688 IF (trim(particles(cur_indices(ip))%atomic_kind%element_symbol) /= trim(elem_mapping(1, ielem))) cycle
697 IF (new_element_extended)
THEN
698 i_rep_section_extended = i_rep_section_extended + 1
700 CALL section_vals_val_set(qm_kind_section,
"_SECTION_PARAMETERS_", i_rep_section=i_rep_section_extended, &
701 c_val=elem_mapping(2, ielem))
702 i_rep_val_extended = 0
703 new_element_extended = .false.
705 i_rep_val_extended = i_rep_val_extended + 1
707 i_rep_val=i_rep_val_extended, i_val=cur_indices(ip))
715 IF (new_element_core)
THEN
716 i_rep_section_core = i_rep_section_core + 1
718 CALL section_vals_val_set(qm_kind_section,
"_SECTION_PARAMETERS_", i_rep_section=i_rep_section_core, &
719 c_val=elem_mapping(2, ielem))
721 new_element_core = .false.
723 i_rep_val_core = i_rep_val_core + 1
725 i_rep_val=i_rep_val_core, i_val=cur_indices(ip))
735 DEALLOCATE (elem_mapping, cur_indices, cur_labels)
737 IF (debug_this_module .AND. output_unit > 0)
THEN
738 WRITE (output_unit, *)
"qmmm_core_section"
740 WRITE (output_unit, *)
"qmmm_extended_section"
752 SUBROUTINE get_force_mixing_indices(force_mixing_section, indices, labels)
753 TYPE(section_vals_type),
POINTER :: force_mixing_section
754 INTEGER,
POINTER :: indices(:), labels(:)
756 INTEGER :: i_rep_val, n_indices, n_labels, n_reps
757 INTEGER,
POINTER :: indices_entry(:), labels_entry(:)
759 TYPE(section_vals_type),
POINTER :: restart_section
761 NULLIFY (indices, labels)
764 IF (.NOT. explicit)
THEN
765 ALLOCATE (indices(0))
773 DO i_rep_val = 1, n_reps
775 i_rep_val=i_rep_val, i_vals=indices_entry)
776 n_indices = n_indices +
SIZE(indices_entry)
778 ALLOCATE (indices(n_indices))
780 DO i_rep_val = 1, n_reps
782 i_rep_val=i_rep_val, i_vals=indices_entry)
783 indices(n_indices + 1:n_indices +
SIZE(indices_entry)) = indices_entry
784 n_indices = n_indices +
SIZE(indices_entry)
789 DO i_rep_val = 1, n_reps
791 i_rep_val=i_rep_val, i_vals=labels_entry)
792 n_labels = n_labels +
SIZE(labels_entry)
794 ALLOCATE (labels(n_labels))
796 DO i_rep_val = 1, n_reps
798 i_rep_val=i_rep_val, i_vals=labels_entry)
799 labels(n_labels + 1:n_labels +
SIZE(labels_entry)) = labels_entry
800 n_labels = n_labels +
SIZE(labels_entry)
803 IF (n_indices /= n_labels) &
804 cpabort(
"got unequal numbers of force_mixing indices and labels!")
805 END SUBROUTINE get_force_mixing_indices
subroutine pbc(r, r_pbc, s, s_pbc, a, b, c, alpha, beta, gamma, debug, info, pbc0, h, hinv)
...
Define the atomic kind types and their sub types.
subroutine, public set_atomic_kind(atomic_kind, element_symbol, name, mass, kind_number, natom, atom_list, fist_potential, shell, shell_active, damping)
Set the components of an atomic kind data set.
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.
Handles all functions related to the CELL.
subroutine, public cell_copy(cell_in, cell_out, tag)
Copy cell variable.
various routines to log and control the output. The idea is that decisions about where to log should ...
recursive integer function, public cp_logger_get_default_unit_nr(logger, local, skip_not_ionode)
asks the default unit number of the given logger. try to use cp_logger_get_unit_nr
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...
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
Define the neighbor list data types and the corresponding functionality.
subroutine, public fist_neighbor_deallocate(fist_neighbor)
...
Generate the atomic neighbor lists for FIST.
subroutine, public build_fist_neighbor_lists(atomic_kind_set, particle_set, local_particles, cell, r_max, r_minsq, ei_scale14, vdw_scale14, nonbonded, para_env, build_from_scratch, geo_check, mm_section, full_nl, exclusions)
...
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
Utility routines for the memory handling.
Interface to the message passing library MPI.
represent a simple array based list of the given type
Define the data structure for the molecule information.
represent a simple array based list of the given type
Define the data structure for the particle information.
logical function, public qmmm_ff_precond_only_qm(id1, id2, id3, id4, is_link)
This function handles the atom names and modifies the "_QM_" prefix, in order to find the parameters ...
integer, parameter, public force_mixing_label_qm_core
integer, parameter, public force_mixing_label_buffer
integer, parameter, public force_mixing_label_none
integer, parameter, public force_mixing_label_qm_core_list
integer, parameter, public force_mixing_label_qm_dynamics_list
integer, parameter, public force_mixing_label_termination
integer, parameter, public force_mixing_label_buffer_list
integer, parameter, public force_mixing_label_qm_dynamics
Basic container type for QM/MM.
subroutine, public qmmm_env_get(qmmm_env, subsys, potential_energy, kinetic_energy)
...
subroutine, public apply_qmmm_translate(qmmm_env)
Apply translation to the full system in order to center the QM system into the QM box.
Basic container type for QM/MM with force mixing.
Routines used for force-mixing QM/MM calculations.
subroutine, public apply_qmmmx_translate(qmmmx_env)
Apply translation to the full system in order to center the QM system into the QM box.
subroutine, public setup_force_mixing_qmmm_sections(subsys, qmmm_section, qmmm_core_section, qmmm_extended_section)
...
subroutine, public update_force_mixing_labels(subsys, qmmm_section, labels_changed)
...