75#include "./base/base_uses.f90"
79 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'force_fields_util'
108 molecule_set, ewald_env, fist_nonbond_env, ff_type, root_section, qmmm, &
109 qmmm_env, mm_section, subsys_section, shell_particle_set, core_particle_set, &
120 LOGICAL,
INTENT(IN),
OPTIONAL :: qmmm
123 TYPE(
particle_type),
DIMENSION(:),
POINTER :: shell_particle_set, core_particle_set
126 CHARACTER(len=*),
PARAMETER :: routinen =
'force_field_pack'
128 CHARACTER(LEN=default_string_length), &
129 DIMENSION(:),
POINTER :: ainfo
130 INTEGER :: handle, iw, iw2, iw3, iw4, output_unit
131 LOGICAL :: do_zbl, explicit, fatal, ignore_fatal, &
133 REAL(kind=
dp) :: ewald_rcut, verlet_skin
134 REAL(kind=
dp),
DIMENSION(:),
POINTER :: charges
143 CALL timeset(routinen, handle)
145 ignore_fatal = ff_type%ignore_missing_critical
146 NULLIFY (logger, ainfo, charges_section, charges)
159 NULLIFY (potparm_nonbond14, potparm_nonbond)
161 IF (
PRESENT(qmmm) .AND.
PRESENT(qmmm_env)) my_qmmm = qmmm
162 inp_info => ff_type%inp_info
163 chm_info => ff_type%chm_info
164 gro_info => ff_type%gro_info
165 amb_info => ff_type%amb_info
199 fatal, ainfo, chm_info, inp_info, gro_info, &
205 fatal, ainfo, chm_info, inp_info, gro_info, &
208 CALL release_ff_missing_par(fatal, ignore_fatal, ainfo, output_unit, iw)
213 ainfo, chm_info, inp_info, iw)
218 ainfo, chm_info, inp_info, gro_info, amb_info, iw)
223 ainfo, chm_info, inp_info, gro_info, iw)
231 CALL release_ff_missing_par(fatal, ignore_fatal, ainfo, output_unit, iw)
235 IF (.NOT. explicit)
THEN
241 ainfo, my_qmmm, inp_info)
243 CALL release_ff_missing_par(fatal, ignore_fatal, ainfo, output_unit, iw)
250 qmmm_env, inp_info, iw4)
265 molecule_kind_set, molecule_set, root_section, subsys_section, &
266 shell_particle_set, core_particle_set, cell, iw, inp_info)
267 IF (ff_type%do_nonbonded)
THEN
273 chm_info, inp_info, gro_info, amb_info, potparm_nonbond14, ewald_env)
275 CALL release_ff_missing_par(fatal, ignore_fatal, ainfo, output_unit, iw)
279 potparm_nonbond14, do_zbl, nonbonded_type=
"NONBONDED14")
285 fatal, iw, ainfo, chm_info, inp_info, gro_info, amb_info, &
286 potparm_nonbond, ewald_env)
288 CALL release_ff_missing_par(fatal, ignore_fatal, ainfo, output_unit, iw)
292 potparm_nonbond, do_zbl, nonbonded_type=
"NONBONDED")
300 ALLOCATE (fist_nonbond_env)
302 potparm_nonbond14, potparm_nonbond, ff_type%do_nonbonded, &
303 ff_type%do_electrostatics, verlet_skin, ewald_rcut, ff_type%ei_scale14, &
304 ff_type%vdw_scale14, ff_type%shift_cutoff)
315 CALL timestop(handle)
327 SUBROUTINE release_ff_missing_par(fatal, ignore_fatal, array, output_unit, iw)
328 LOGICAL,
INTENT(INOUT) :: fatal, ignore_fatal
329 CHARACTER(LEN=default_string_length), &
330 DIMENSION(:),
POINTER :: array
331 INTEGER,
INTENT(IN) :: output_unit, iw
335 IF (
ASSOCIATED(array))
THEN
336 IF (output_unit > 0)
THEN
337 WRITE (output_unit, *)
338 WRITE (output_unit,
'(T2,"FORCEFIELD|",A)') &
339 " WARNING: A non Critical ForceField parameter is missing! CP2K GOES!", &
340 " Non critical parameters are:Urey-Bradley,Torsions,Impropers, Opbends and 1-4", &
341 " All missing parameters will not contribute to the potential energy!"
342 IF (fatal .OR. iw > 0)
THEN
343 WRITE (output_unit, *)
344 DO i = 1,
SIZE(array)
345 WRITE (output_unit,
'(A)') array(i)
348 IF (.NOT. fatal .AND. iw < 0)
THEN
349 WRITE (output_unit,
'(T2,"FORCEFIELD|",A)') &
350 " Activate the print key FF_INFO to have a list of missing parameters"
351 WRITE (output_unit, *)
357 IF (ignore_fatal)
THEN
358 IF (output_unit > 0)
THEN
359 WRITE (output_unit, *)
360 WRITE (output_unit,
'(T2,"FORCEFIELD|",A)') &
361 " WARNING: Ignoring missing critical FF parameters! CP2K GOES!", &
362 " Critical parameters are: Bonds, Bends and Charges", &
363 " All missing parameters will not contribute to the potential energy!"
366 cpabort(
"Missing critical ForceField parameters!")
369 END SUBROUTINE release_ff_missing_par
380 molecule_set, mm_section, charges)
386 REAL(kind=
dp),
DIMENSION(:),
POINTER :: charges
388 CHARACTER(len=*),
PARAMETER :: routinen =
'force_field_qeff_output'
390 CHARACTER(LEN=default_string_length) :: atmname, molname
391 INTEGER :: first, handle, iatom, imol, iw, j, jatom
392 LOGICAL :: shell_active
393 REAL(kind=
dp) :: mass, mass_mol, mass_sum, qeff, &
395 TYPE(
atom_type),
DIMENSION(:),
POINTER :: atom_list
402 CALL timeset(routinen, handle)
416 DO imol = 1,
SIZE(molecule_kind_set)
419 molecule_kind => molecule_kind_set(imol)
421 j = molecule_kind_set(imol)%molecule_list(1)
422 molecule => molecule_set(j)
426 name=molname, atom_list=atom_list)
427 DO iatom = 1,
SIZE(atom_list)
428 atomic_kind => atom_list(iatom)%atomic_kind
430 name=atmname, qeff=qeff, mass=mass, shell_active=shell_active, shell=shell)
431 IF (shell_active) qeff = shell%charge_core + shell%charge_shell
432 IF (
ASSOCIATED(charges))
THEN
433 jatom = first - 1 + iatom
434 qeff = charges(jatom)
436 IF (iw > 0)
WRITE (iw, *)
" atom ", iatom,
" ", trim(atmname),
" charge = ", qeff
437 qeff_mol = qeff_mol + qeff
438 mass_mol = mass_mol + mass
440 CALL set_molecule_kind(molecule_kind=molecule_kind, charge=qeff_mol, mass=mass_mol)
441 IF (iw > 0)
WRITE (iw, *)
" Mol Kind ", trim(molname),
" charge = ", qeff_mol
446 DO iatom = 1,
SIZE(particle_set)
447 atomic_kind => particle_set(iatom)%atomic_kind
449 name=atmname, qeff=qeff, mass=mass, shell_active=shell_active, shell=shell)
450 IF (shell_active) qeff = shell%charge_core + shell%charge_shell
451 IF (
ASSOCIATED(charges))
THEN
452 qeff = charges(iatom)
454 IF (iw > 0)
WRITE (iw, *)
" atom ", iatom,
" ", trim(atmname), &
456 qeff_sum = qeff_sum + qeff
457 mass_sum = mass_sum + mass
459 IF (iw > 0)
WRITE (iw,
'(A,F20.10)')
" Total system charge = ", qeff_sum
460 IF (iw > 0)
WRITE (iw,
'(A,F20.10)')
" Total system mass = ", mass_sum
465 CALL timestop(handle)
479 CHARACTER(len=*),
PARAMETER :: routinen =
'clean_intra_force_kind'
481 INTEGER :: atm2_a, atm2_b, atm2_c, atm_a, atm_b, atm_c, atm_d, counter, handle, i, ibend, &
482 ibond, icolv, ig3x3, ig4x6, iimpr, ikind, iopbend, itorsion, iub, iw, j, k, nbend, nbond, &
483 newkind, ng3x3, ng4x6, nimpr, nopbend, ntorsion, nub
484 INTEGER,
POINTER :: bad1(:), bad2(:)
485 LOGICAL :: unsetme, valid_kind
486 TYPE(
bend_kind_type),
DIMENSION(:),
POINTER :: bend_kind_set, new_bend_kind_set
487 TYPE(
bend_type),
DIMENSION(:),
POINTER :: bend_list, new_bend_list
488 TYPE(
bond_kind_type),
DIMENSION(:),
POINTER :: bond_kind_set, new_bond_kind_set
489 TYPE(
bond_type),
DIMENSION(:),
POINTER :: bond_list, new_bond_list
495 TYPE(
impr_kind_type),
DIMENSION(:),
POINTER :: impr_kind_set, new_impr_kind_set
496 TYPE(
impr_type),
DIMENSION(:),
POINTER :: impr_list, new_impr_list
498 TYPE(
opbend_kind_type),
DIMENSION(:),
POINTER :: new_opbend_kind_set, opbend_kind_set
499 TYPE(
opbend_type),
DIMENSION(:),
POINTER :: new_opbend_list, opbend_list
500 TYPE(
torsion_kind_type),
DIMENSION(:),
POINTER :: new_torsion_kind_set, torsion_kind_set
501 TYPE(
torsion_type),
DIMENSION(:),
POINTER :: new_torsion_list, torsion_list
502 TYPE(
ub_kind_type),
DIMENSION(:),
POINTER :: new_ub_kind_set, ub_kind_set
503 TYPE(
ub_type),
DIMENSION(:),
POINTER :: new_ub_list, ub_list
505 CALL timeset(routinen, handle)
513 DO ikind = 1,
SIZE(molecule_kind_set)
514 molecule_kind => molecule_kind_set(ikind)
516 colv_list=colv_list, &
519 IF (
ASSOCIATED(colv_list))
THEN
520 DO icolv = 1,
SIZE(colv_list)
522 ((.NOT. colv_list(icolv)%use_points) .OR. (
SIZE(colv_list(icolv)%i_atoms) == 2)))
THEN
523 atm_a = colv_list(icolv)%i_atoms(1)
524 atm_b = colv_list(icolv)%i_atoms(2)
527 atm2_a = bond_list(ibond)%a
528 atm2_b = bond_list(ibond)%b
529 IF (atm2_a == atm_a .AND. atm2_b == atm_b) unsetme = .true.
530 IF (atm2_a == atm_b .AND. atm2_b == atm_a) unsetme = .true.
531 IF (unsetme) bond_list(ibond)%id_type =
do_ff_undef
540 DO ikind = 1,
SIZE(molecule_kind_set)
541 molecule_kind => molecule_kind_set(ikind)
543 colv_list=colv_list, &
546 IF (
ASSOCIATED(colv_list))
THEN
549 atm_a = bend_list(ibend)%a
550 atm_b = bend_list(ibend)%b
551 atm_c = bend_list(ibend)%c
552 DO icolv = 1,
SIZE(colv_list)
554 ((.NOT. colv_list(icolv)%use_points) .OR. (
SIZE(colv_list(icolv)%i_atoms) == 2)))
THEN
555 atm2_a = colv_list(icolv)%i_atoms(1)
556 atm2_b = colv_list(icolv)%i_atoms(2)
559 IF (((atm2_a == atm_a) .AND. (atm2_b == atm_c)) .OR. &
560 ((atm2_a == atm_c) .AND. (atm2_b == atm_a)))
THEN
566 IF (unsetme) bend_list(ibend)%id_type =
do_ff_undef
573 DO ikind = 1,
SIZE(molecule_kind_set)
574 molecule_kind => molecule_kind_set(ikind)
577 g3x3_list=g3x3_list, &
581 atm_a = g3x3_list(ig3x3)%a
582 atm_b = g3x3_list(ig3x3)%b
583 atm_c = g3x3_list(ig3x3)%c
586 atm2_a = bond_list(ibond)%a
587 atm2_b = bond_list(ibond)%b
588 IF (atm2_a == atm_a .AND. atm2_b == atm_b) unsetme = .true.
589 IF (atm2_a == atm_a .AND. atm2_b == atm_c) unsetme = .true.
590 IF (atm2_a == atm_c .AND. atm2_b == atm_c) unsetme = .true.
591 IF (unsetme) bond_list(ibond)%id_type =
do_ff_undef
598 DO ikind = 1,
SIZE(molecule_kind_set)
599 molecule_kind => molecule_kind_set(ikind)
602 g3x3_list=g3x3_list, &
606 atm_a = g3x3_list(ig3x3)%a
607 atm_b = g3x3_list(ig3x3)%b
608 atm_c = g3x3_list(ig3x3)%c
611 atm2_a = bend_list(ibend)%a
612 atm2_b = bend_list(ibend)%b
613 atm2_c = bend_list(ibend)%c
614 IF (atm2_a == atm_a .AND. atm2_b == atm_b .AND. atm2_c == atm_c) unsetme = .true.
615 IF (atm2_a == atm_a .AND. atm2_b == atm_c .AND. atm2_c == atm_b) unsetme = .true.
616 IF (atm2_a == atm_b .AND. atm2_b == atm_a .AND. atm2_c == atm_c) unsetme = .true.
617 IF (atm2_a == atm_b .AND. atm2_b == atm_c .AND. atm2_c == atm_a) unsetme = .true.
618 IF (unsetme) bend_list(ibend)%id_type =
do_ff_undef
625 DO ikind = 1,
SIZE(molecule_kind_set)
626 molecule_kind => molecule_kind_set(ikind)
629 g4x6_list=g4x6_list, &
633 atm_a = g4x6_list(ig4x6)%a
634 atm_b = g4x6_list(ig4x6)%b
635 atm_c = g4x6_list(ig4x6)%c
636 atm_d = g4x6_list(ig4x6)%d
639 atm2_a = bond_list(ibond)%a
640 atm2_b = bond_list(ibond)%b
641 IF (atm2_a == atm_a .AND. atm2_b == atm_b) unsetme = .true.
642 IF (atm2_a == atm_a .AND. atm2_b == atm_c) unsetme = .true.
643 IF (atm2_a == atm_a .AND. atm2_b == atm_d) unsetme = .true.
644 IF (unsetme) bond_list(ibond)%id_type =
do_ff_undef
651 DO ikind = 1,
SIZE(molecule_kind_set)
652 molecule_kind => molecule_kind_set(ikind)
655 g4x6_list=g4x6_list, &
659 atm_a = g4x6_list(ig4x6)%a
660 atm_b = g4x6_list(ig4x6)%b
661 atm_c = g4x6_list(ig4x6)%c
662 atm_d = g4x6_list(ig4x6)%d
665 atm2_a = bend_list(ibend)%a
666 atm2_b = bend_list(ibend)%b
667 atm2_c = bend_list(ibend)%c
668 IF (atm2_a == atm_b .AND. atm2_b == atm_a .AND. atm2_c == atm_c) unsetme = .true.
669 IF (atm2_a == atm_b .AND. atm2_b == atm_a .AND. atm2_c == atm_d) unsetme = .true.
670 IF (atm2_a == atm_c .AND. atm2_b == atm_a .AND. atm2_c == atm_d) unsetme = .true.
671 IF (unsetme) bend_list(ibend)%id_type =
do_ff_undef
678 DO ikind = 1,
SIZE(molecule_kind_set)
679 molecule_kind => molecule_kind_set(ikind)
682 bond_kind_set=bond_kind_set, &
685 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") Old BOND Count: ", &
686 SIZE(bond_list),
SIZE(bond_kind_set)
687 IF (iw > 0)
WRITE (iw,
'(2I6)') (bond_list(ibond)%a, bond_list(ibond)%b, ibond=1,
SIZE(bond_list))
689 ALLOCATE (bad1(
SIZE(bond_kind_set)))
691 DO ibond = 1,
SIZE(bond_kind_set)
693 IF (bond_kind_set(ibond)%id_type ==
do_ff_undef) unsetme = .true.
695 DO i = 1,
SIZE(bond_list)
697 bond_list(i)%bond_kind%kind_number == ibond)
THEN
702 IF (.NOT. valid_kind) unsetme = .true.
703 IF (unsetme) bad1(ibond) = 1
705 IF (sum(bad1) /= 0)
THEN
706 counter =
SIZE(bond_kind_set) - sum(bad1)
709 DO ibond = 1,
SIZE(bond_kind_set)
710 IF (bad1(ibond) == 0)
THEN
711 counter = counter + 1
712 new_bond_kind_set(counter) = bond_kind_set(ibond)
716 ALLOCATE (bad2(
SIZE(bond_list)))
718 DO ibond = 1,
SIZE(bond_list)
720 IF (bond_list(ibond)%bond_kind%id_type ==
do_ff_undef) unsetme = .true.
721 IF (bond_list(ibond)%id_type ==
do_ff_undef) unsetme = .true.
722 IF (unsetme) bad2(ibond) = 1
724 IF (sum(bad2) /= 0)
THEN
725 counter =
SIZE(bond_list) - sum(bad2)
726 ALLOCATE (new_bond_list(counter))
728 DO ibond = 1,
SIZE(bond_list)
729 IF (bad2(ibond) == 0)
THEN
730 counter = counter + 1
731 new_bond_list(counter) = bond_list(ibond)
732 newkind = bond_list(ibond)%bond_kind%kind_number
733 newkind = newkind - sum(bad1(1:newkind))
734 new_bond_list(counter)%bond_kind => new_bond_kind_set(newkind)
738 nbond=
SIZE(new_bond_list), &
739 bond_kind_set=new_bond_kind_set, &
740 bond_list=new_bond_list)
741 DO ibond = 1,
SIZE(new_bond_kind_set)
742 new_bond_kind_set(ibond)%kind_number = ibond
747 DEALLOCATE (bond_list)
748 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") New BOND Count: ", &
749 SIZE(new_bond_list),
SIZE(new_bond_kind_set)
750 IF (iw > 0)
WRITE (iw,
'(2I6)') (new_bond_list(ibond)%a, new_bond_list(ibond)%b, &
751 ibond=1,
SIZE(new_bond_list))
759 DO ikind = 1,
SIZE(molecule_kind_set)
760 molecule_kind => molecule_kind_set(ikind)
763 bend_kind_set=bend_kind_set, &
766 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") Old BEND Count: ", &
767 SIZE(bend_list),
SIZE(bend_kind_set)
768 IF (iw > 0)
WRITE (iw,
'(3I6)') (bend_list(ibend)%a, bend_list(ibend)%b, &
769 bend_list(ibend)%c, ibend=1,
SIZE(bend_list))
771 ALLOCATE (bad1(
SIZE(bend_kind_set)))
773 DO ibend = 1,
SIZE(bend_kind_set)
775 IF (bend_kind_set(ibend)%id_type ==
do_ff_undef) unsetme = .true.
777 DO i = 1,
SIZE(bend_list)
779 bend_list(i)%bend_kind%kind_number == ibend)
THEN
784 IF (.NOT. valid_kind) unsetme = .true.
785 IF (unsetme) bad1(ibend) = 1
787 IF (sum(bad1) /= 0)
THEN
788 counter =
SIZE(bend_kind_set) - sum(bad1)
791 DO ibend = 1,
SIZE(bend_kind_set)
792 IF (bad1(ibend) == 0)
THEN
793 counter = counter + 1
794 new_bend_kind_set(counter) = bend_kind_set(ibend)
798 ALLOCATE (bad2(
SIZE(bend_list)))
800 DO ibend = 1,
SIZE(bend_list)
802 IF (bend_list(ibend)%bend_kind%id_type ==
do_ff_undef) unsetme = .true.
803 IF (bend_list(ibend)%id_type ==
do_ff_undef) unsetme = .true.
804 IF (unsetme) bad2(ibend) = 1
806 IF (sum(bad2) /= 0)
THEN
807 counter =
SIZE(bend_list) - sum(bad2)
808 ALLOCATE (new_bend_list(counter))
810 DO ibend = 1,
SIZE(bend_list)
811 IF (bad2(ibend) == 0)
THEN
812 counter = counter + 1
813 new_bend_list(counter) = bend_list(ibend)
814 newkind = bend_list(ibend)%bend_kind%kind_number
815 newkind = newkind - sum(bad1(1:newkind))
816 new_bend_list(counter)%bend_kind => new_bend_kind_set(newkind)
820 nbend=
SIZE(new_bend_list), &
821 bend_kind_set=new_bend_kind_set, &
822 bend_list=new_bend_list)
823 DO ibend = 1,
SIZE(new_bend_kind_set)
824 new_bend_kind_set(ibend)%kind_number = ibend
829 DEALLOCATE (bend_list)
830 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") New BEND Count: ", &
831 SIZE(new_bend_list),
SIZE(new_bend_kind_set)
832 IF (iw > 0)
WRITE (iw,
'(3I6)') (new_bend_list(ibend)%a, new_bend_list(ibend)%b, &
833 new_bend_list(ibend)%c, ibend=1,
SIZE(new_bend_list))
841 DO ikind = 1,
SIZE(molecule_kind_set)
842 molecule_kind => molecule_kind_set(ikind)
845 ub_kind_set=ub_kind_set, &
848 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") Old UB Count: ", &
849 SIZE(ub_list),
SIZE(ub_kind_set)
850 IF (iw > 0)
WRITE (iw,
'(3I6)') (ub_list(iub)%a, ub_list(iub)%b, &
851 ub_list(iub)%c, iub=1,
SIZE(ub_list))
853 ALLOCATE (bad1(
SIZE(ub_kind_set)))
855 DO iub = 1,
SIZE(ub_kind_set)
857 IF (ub_kind_set(iub)%id_type ==
do_ff_undef) unsetme = .true.
859 DO i = 1,
SIZE(ub_list)
861 ub_list(i)%ub_kind%kind_number == iub)
THEN
866 IF (.NOT. valid_kind) unsetme = .true.
867 IF (unsetme) bad1(iub) = 1
869 IF (sum(bad1) /= 0)
THEN
870 counter =
SIZE(ub_kind_set) - sum(bad1)
873 DO iub = 1,
SIZE(ub_kind_set)
874 IF (bad1(iub) == 0)
THEN
875 counter = counter + 1
876 new_ub_kind_set(counter) = ub_kind_set(iub)
880 ALLOCATE (bad2(
SIZE(ub_list)))
882 DO iub = 1,
SIZE(ub_list)
884 IF (ub_list(iub)%ub_kind%id_type ==
do_ff_undef) unsetme = .true.
885 IF (ub_list(iub)%id_type ==
do_ff_undef) unsetme = .true.
886 IF (unsetme) bad2(iub) = 1
888 IF (sum(bad2) /= 0)
THEN
889 counter =
SIZE(ub_list) - sum(bad2)
890 ALLOCATE (new_ub_list(counter))
892 DO iub = 1,
SIZE(ub_list)
893 IF (bad2(iub) == 0)
THEN
894 counter = counter + 1
895 new_ub_list(counter) = ub_list(iub)
896 newkind = ub_list(iub)%ub_kind%kind_number
897 newkind = newkind - sum(bad1(1:newkind))
898 new_ub_list(counter)%ub_kind => new_ub_kind_set(newkind)
902 nub=
SIZE(new_ub_list), &
903 ub_kind_set=new_ub_kind_set, &
905 DO iub = 1,
SIZE(new_ub_kind_set)
906 new_ub_kind_set(iub)%kind_number = iub
912 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") New UB Count: ", &
913 SIZE(new_ub_list),
SIZE(new_ub_kind_set)
914 IF (iw > 0)
WRITE (iw,
'(3I6)') (new_ub_list(iub)%a, new_ub_list(iub)%b, &
915 new_ub_list(iub)%c, iub=1,
SIZE(new_ub_list))
923 DO ikind = 1,
SIZE(molecule_kind_set)
924 molecule_kind => molecule_kind_set(ikind)
927 torsion_kind_set=torsion_kind_set, &
928 torsion_list=torsion_list)
929 IF (ntorsion > 0)
THEN
930 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") Old TORSION Count: ", &
931 SIZE(torsion_list),
SIZE(torsion_kind_set)
932 IF (iw > 0)
WRITE (iw,
'(4I6)') (torsion_list(itorsion)%a, torsion_list(itorsion)%b, &
933 torsion_list(itorsion)%c, torsion_list(itorsion)%d, itorsion=1,
SIZE(torsion_list))
935 ALLOCATE (bad1(
SIZE(torsion_kind_set)))
937 DO itorsion = 1,
SIZE(torsion_kind_set)
939 IF (torsion_kind_set(itorsion)%id_type ==
do_ff_undef) unsetme = .true.
941 DO i = 1,
SIZE(torsion_list)
943 torsion_list(i)%torsion_kind%kind_number == itorsion)
THEN
948 IF (.NOT. valid_kind) unsetme = .true.
949 IF (unsetme) bad1(itorsion) = 1
951 IF (sum(bad1) /= 0)
THEN
952 counter =
SIZE(torsion_kind_set) - sum(bad1)
955 DO itorsion = 1,
SIZE(torsion_kind_set)
956 IF (bad1(itorsion) == 0)
THEN
957 counter = counter + 1
958 new_torsion_kind_set(counter) = torsion_kind_set(itorsion)
959 i =
SIZE(torsion_kind_set(itorsion)%m)
960 j =
SIZE(torsion_kind_set(itorsion)%k)
961 k =
SIZE(torsion_kind_set(itorsion)%phi0)
962 ALLOCATE (new_torsion_kind_set(counter)%m(i))
963 ALLOCATE (new_torsion_kind_set(counter)%k(i))
964 ALLOCATE (new_torsion_kind_set(counter)%phi0(i))
965 new_torsion_kind_set(counter)%m = torsion_kind_set(itorsion)%m
966 new_torsion_kind_set(counter)%k = torsion_kind_set(itorsion)%k
967 new_torsion_kind_set(counter)%phi0 = torsion_kind_set(itorsion)%phi0
971 ALLOCATE (bad2(
SIZE(torsion_list)))
973 DO itorsion = 1,
SIZE(torsion_list)
975 IF (torsion_list(itorsion)%torsion_kind%id_type ==
do_ff_undef) unsetme = .true.
976 IF (torsion_list(itorsion)%id_type ==
do_ff_undef) unsetme = .true.
977 IF (unsetme) bad2(itorsion) = 1
979 IF (sum(bad2) /= 0)
THEN
980 counter =
SIZE(torsion_list) - sum(bad2)
981 ALLOCATE (new_torsion_list(counter))
983 DO itorsion = 1,
SIZE(torsion_list)
984 IF (bad2(itorsion) == 0)
THEN
985 counter = counter + 1
986 new_torsion_list(counter) = torsion_list(itorsion)
987 newkind = torsion_list(itorsion)%torsion_kind%kind_number
988 newkind = newkind - sum(bad1(1:newkind))
989 new_torsion_list(counter)%torsion_kind => new_torsion_kind_set(newkind)
993 ntorsion=
SIZE(new_torsion_list), &
994 torsion_kind_set=new_torsion_kind_set, &
995 torsion_list=new_torsion_list)
996 DO itorsion = 1,
SIZE(new_torsion_kind_set)
997 new_torsion_kind_set(itorsion)%kind_number = itorsion
1001 DO itorsion = 1,
SIZE(torsion_kind_set)
1004 DEALLOCATE (torsion_kind_set)
1005 DEALLOCATE (torsion_list)
1006 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") New TORSION Count: ", &
1007 SIZE(new_torsion_list),
SIZE(new_torsion_kind_set)
1008 IF (iw > 0)
WRITE (iw,
'(4I6)') (new_torsion_list(itorsion)%a, new_torsion_list(itorsion)%b, &
1009 new_torsion_list(itorsion)%c, new_torsion_list(itorsion)%d, itorsion=1, &
1010 SIZE(new_torsion_list))
1018 DO ikind = 1,
SIZE(molecule_kind_set)
1019 molecule_kind => molecule_kind_set(ikind)
1022 impr_kind_set=impr_kind_set, &
1023 impr_list=impr_list)
1025 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") Old IMPROPER Count: ", &
1026 SIZE(impr_list),
SIZE(impr_kind_set)
1027 NULLIFY (bad1, bad2)
1028 ALLOCATE (bad1(
SIZE(impr_kind_set)))
1030 DO iimpr = 1,
SIZE(impr_kind_set)
1032 IF (impr_kind_set(iimpr)%id_type ==
do_ff_undef) unsetme = .true.
1033 valid_kind = .false.
1034 DO i = 1,
SIZE(impr_list)
1036 impr_list(i)%impr_kind%kind_number == iimpr)
THEN
1041 IF (.NOT. valid_kind) unsetme = .true.
1042 IF (unsetme) bad1(iimpr) = 1
1044 IF (sum(bad1) /= 0)
THEN
1045 counter =
SIZE(impr_kind_set) - sum(bad1)
1048 DO iimpr = 1,
SIZE(impr_kind_set)
1049 IF (bad1(iimpr) == 0)
THEN
1050 counter = counter + 1
1051 new_impr_kind_set(counter) = impr_kind_set(iimpr)
1055 ALLOCATE (bad2(
SIZE(impr_list)))
1057 DO iimpr = 1,
SIZE(impr_list)
1059 IF (impr_list(iimpr)%impr_kind%id_type ==
do_ff_undef) unsetme = .true.
1060 IF (impr_list(iimpr)%id_type ==
do_ff_undef) unsetme = .true.
1061 IF (unsetme) bad2(iimpr) = 1
1063 IF (sum(bad2) /= 0)
THEN
1064 counter =
SIZE(impr_list) - sum(bad2)
1065 ALLOCATE (new_impr_list(counter))
1067 DO iimpr = 1,
SIZE(impr_list)
1068 IF (bad2(iimpr) == 0)
THEN
1069 counter = counter + 1
1070 new_impr_list(counter) = impr_list(iimpr)
1071 newkind = impr_list(iimpr)%impr_kind%kind_number
1072 newkind = newkind - sum(bad1(1:newkind))
1073 new_impr_list(counter)%impr_kind => new_impr_kind_set(newkind)
1077 nimpr=
SIZE(new_impr_list), &
1078 impr_kind_set=new_impr_kind_set, &
1079 impr_list=new_impr_list)
1080 DO iimpr = 1,
SIZE(new_impr_kind_set)
1081 new_impr_kind_set(iimpr)%kind_number = iimpr
1085 DO iimpr = 1,
SIZE(impr_kind_set)
1088 DEALLOCATE (impr_kind_set)
1089 DEALLOCATE (impr_list)
1090 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") New IMPROPER Count: ", &
1091 SIZE(new_impr_list),
SIZE(new_impr_kind_set)
1099 DO ikind = 1,
SIZE(molecule_kind_set)
1100 molecule_kind => molecule_kind_set(ikind)
1103 opbend_kind_set=opbend_kind_set, &
1104 opbend_list=opbend_list)
1105 IF (nopbend > 0)
THEN
1106 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") Old OPBEND Count: ", &
1107 SIZE(opbend_list),
SIZE(opbend_kind_set)
1108 NULLIFY (bad1, bad2)
1109 ALLOCATE (bad1(
SIZE(opbend_kind_set)))
1111 DO iopbend = 1,
SIZE(opbend_kind_set)
1113 IF (opbend_kind_set(iopbend)%id_type ==
do_ff_undef) unsetme = .true.
1114 valid_kind = .false.
1115 DO i = 1,
SIZE(opbend_list)
1117 opbend_list(i)%opbend_kind%kind_number == iopbend)
THEN
1122 IF (.NOT. valid_kind) unsetme = .true.
1123 IF (unsetme) bad1(iopbend) = 1
1125 IF (sum(bad1) /= 0)
THEN
1126 counter =
SIZE(opbend_kind_set) - sum(bad1)
1129 DO iopbend = 1,
SIZE(opbend_kind_set)
1130 IF (bad1(iopbend) == 0)
THEN
1131 counter = counter + 1
1132 new_opbend_kind_set(counter) = opbend_kind_set(iopbend)
1136 ALLOCATE (bad2(
SIZE(opbend_list)))
1138 DO iopbend = 1,
SIZE(opbend_list)
1140 IF (opbend_list(iopbend)%opbend_kind%id_type ==
do_ff_undef) unsetme = .true.
1141 IF (opbend_list(iopbend)%id_type ==
do_ff_undef) unsetme = .true.
1142 IF (unsetme) bad2(iopbend) = 1
1144 IF (sum(bad2) /= 0)
THEN
1145 counter =
SIZE(opbend_list) - sum(bad2)
1146 ALLOCATE (new_opbend_list(counter))
1148 DO iopbend = 1,
SIZE(opbend_list)
1149 IF (bad2(iopbend) == 0)
THEN
1150 counter = counter + 1
1151 new_opbend_list(counter) = opbend_list(iopbend)
1152 newkind = opbend_list(iopbend)%opbend_kind%kind_number
1153 newkind = newkind - sum(bad1(1:newkind))
1154 new_opbend_list(counter)%opbend_kind => new_opbend_kind_set(newkind)
1158 nopbend=
SIZE(new_opbend_list), &
1159 opbend_kind_set=new_opbend_kind_set, &
1160 opbend_list=new_opbend_list)
1161 DO iopbend = 1,
SIZE(new_opbend_kind_set)
1162 new_opbend_kind_set(iopbend)%kind_number = iopbend
1166 DEALLOCATE (opbend_kind_set)
1167 DEALLOCATE (opbend_list)
1168 IF (iw > 0)
WRITE (iw, *)
" Mol(", ikind,
") New OPBEND Count: ", &
1169 SIZE(new_opbend_list),
SIZE(new_opbend_kind_set)
1181 CALL timestop(handle)
1198 var_values, size_variables, i_rep_sec, input_variables)
1200 CHARACTER(LEN=*),
INTENT(IN) :: func_name
1201 CHARACTER(LEN=default_path_length),
INTENT(OUT) :: xfunction
1202 CHARACTER(LEN=default_string_length), &
1203 DIMENSION(:),
POINTER :: parameters
1204 REAL(kind=
dp),
DIMENSION(:),
POINTER :: values
1205 REAL(kind=
dp),
DIMENSION(:),
OPTIONAL,
POINTER :: var_values
1206 INTEGER,
INTENT(IN),
OPTIONAL :: size_variables, i_rep_sec
1207 CHARACTER(LEN=*),
DIMENSION(:),
OPTIONAL :: input_variables
1209 CHARACTER(LEN=default_string_length), &
1210 DIMENSION(:),
POINTER :: my_par, my_par_tmp, my_units, &
1211 my_units_tmp, my_var
1212 INTEGER :: i, ind, irep, isize, j, mydim, n_par, &
1213 n_units, n_val, nblank
1215 REAL(kind=
dp),
DIMENSION(:),
POINTER :: my_val, my_val_tmp
1217 NULLIFY (my_var, my_par, my_val, my_par_tmp, my_val_tmp)
1219 NULLIFY (my_units_tmp)
1220 IF (
ASSOCIATED(parameters))
THEN
1221 DEALLOCATE (parameters)
1223 IF (
ASSOCIATED(values))
THEN
1227 IF (
PRESENT(i_rep_sec)) irep = i_rep_sec
1230 CALL compress(xfunction, full=.true.)
1231 IF (
PRESENT(input_variables))
THEN
1232 ALLOCATE (my_var(
SIZE(input_variables)))
1233 my_var = input_variables
1237 IF (
ASSOCIATED(my_var))
THEN
1238 mydim =
SIZE(my_var)
1240 IF (
PRESENT(size_variables))
THEN
1241 cpassert(mydim == size_variables)
1246 check = (n_par > 0) .EQV. (n_val > 0)
1251 ALLOCATE (my_par(0))
1252 ALLOCATE (my_val(0))
1254 isize =
SIZE(my_par)
1255 CALL section_vals_val_get(gen_section,
"PARAMETERS", i_rep_section=irep, i_rep_val=i, c_vals=my_par_tmp)
1256 nblank = count(my_par_tmp ==
"")
1257 CALL reallocate(my_par, 1, isize +
SIZE(my_par_tmp) - nblank)
1259 DO j = 1,
SIZE(my_par_tmp)
1260 IF (my_par_tmp(j) ==
"") cycle
1262 my_par(isize + ind) = my_par_tmp(j)
1266 isize =
SIZE(my_val)
1267 CALL section_vals_val_get(gen_section,
"VALUES", i_rep_section=irep, i_rep_val=i, r_vals=my_val_tmp)
1268 CALL reallocate(my_val, 1, isize +
SIZE(my_val_tmp))
1269 my_val(isize + 1:isize +
SIZE(my_val_tmp)) = my_val_tmp
1271 cpassert(
SIZE(my_par) ==
SIZE(my_val))
1273 ALLOCATE (my_units(0))
1274 IF (n_units > 0)
THEN
1276 isize =
SIZE(my_units)
1277 CALL section_vals_val_get(gen_section,
"UNITS", i_rep_section=irep, i_rep_val=i, c_vals=my_units_tmp)
1278 nblank = count(my_units_tmp ==
"")
1279 CALL reallocate(my_units, 1, isize +
SIZE(my_units_tmp) - nblank)
1281 DO j = 1,
SIZE(my_units_tmp)
1282 IF (my_units_tmp(j) ==
"") cycle
1284 my_units(isize + ind) = my_units_tmp(j)
1287 cpassert(
SIZE(my_units) ==
SIZE(my_val))
1289 mydim = mydim +
SIZE(my_val)
1290 IF (
SIZE(my_val) == 0)
THEN
1293 DEALLOCATE (my_units)
1297 ALLOCATE (parameters(mydim))
1298 ALLOCATE (values(mydim))
1300 parameters(1:
SIZE(my_var)) = my_var
1301 values(1:
SIZE(my_var)) = 0.0_dp
1302 IF (
PRESENT(var_values))
THEN
1303 cpassert(
SIZE(var_values) ==
SIZE(my_var))
1304 values(1:
SIZE(my_var)) = var_values
1306 IF (
ASSOCIATED(my_val))
THEN
1307 DO i = 1,
SIZE(my_val)
1308 parameters(
SIZE(my_var) + i) = my_par(i)
1309 IF (n_units > 0)
THEN
1310 values(
SIZE(my_var) + i) =
cp_unit_to_cp2k(my_val(i), trim(adjustl(my_units(i))))
1312 values(
SIZE(my_var) + i) = my_val(i)
1317 IF (
ASSOCIATED(my_par))
THEN
1320 IF (
ASSOCIATED(my_val))
THEN
1323 IF (
ASSOCIATED(my_units))
THEN
1324 DEALLOCATE (my_units)
1326 IF (
PRESENT(input_variables))
THEN
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.
Handles all functions related to the CELL.
Initialize the collective variables types.
integer, parameter, public dist_colvar_id
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
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,...
real(kind=dp) function, public cp_unit_to_cp2k(value, unit_str, defaults, power)
converts to the internal cp2k units to the given unit
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 fist_nonbond_env_create(fist_nonbond_env, atomic_kind_set, potparm14, potparm, do_nonbonded, do_electrostatics, verlet_skin, ewald_rcut, ei_scale14, vdw_scale14, shift_cutoff)
allocates and intitializes a fist_nonbond_env
subroutine, public fist_nonbond_env_set(fist_nonbond_env, potparm14, potparm, rlist_cut, rlist_lowsq, nonbonded, aup, lup, ei_scale14, vdw_scale14, shift_cutoff, do_electrostatics, r_last_update, r_last_update_pbc, rshell_last_update_pbc, rcore_last_update_pbc, cell_last_update, num_update, last_update, counter, natom_types, long_range_correction, eam_data, quip_data, nequip_data, allegro_data, deepmd_data, ace_data, charges)
sets a fist_nonbond_env
Define all structure types related to force field kinds.
pure subroutine, public deallocate_bend_kind_set(bend_kind_set)
Deallocate a bend kind set.
integer, parameter, public do_ff_undef
pure subroutine, public allocate_impr_kind_set(impr_kind_set, nkind)
Allocate and initialize a impr kind set.
pure subroutine, public torsion_kind_dealloc_ref(torsion_kind)
Deallocate a torsion kind element.
pure subroutine, public allocate_torsion_kind_set(torsion_kind_set, nkind)
Allocate and initialize a torsion kind set.
pure subroutine, public allocate_bond_kind_set(bond_kind_set, nkind)
Allocate and initialize a bond kind set.
pure subroutine, public allocate_opbend_kind_set(opbend_kind_set, nkind)
Allocate and initialize a opbend kind set.
pure subroutine, public ub_kind_dealloc_ref(ub_kind_set)
Deallocate a ub kind set.
pure subroutine, public allocate_bend_kind_set(bend_kind_set, nkind)
Allocate and initialize a bend kind set.
pure subroutine, public deallocate_bond_kind_set(bond_kind_set)
Deallocate a bond kind set.
pure subroutine, public allocate_ub_kind_set(ub_kind_set, nkind)
Allocate and initialize a ub kind set.
pure subroutine, public impr_kind_dealloc_ref()
Deallocate a impr kind element.
Define all structures types related to force_fields.
subroutine, public force_field_unique_ub(particle_set, molecule_kind_set, molecule_set, iw)
Determine the number of unique Urey-Bradley kind and allocate ub_kind_set.
subroutine, public force_field_pack_tors(particle_set, molecule_kind_set, molecule_set, ainfo, chm_info, inp_info, gro_info, amb_info, iw)
Pack in torsion information needed for the force_field.
subroutine, public force_field_unique_bond(particle_set, molecule_kind_set, molecule_set, ff_type, iw)
Determine the number of unique bond kind and allocate bond_kind_set.
subroutine, public force_field_unique_opbend(particle_set, molecule_kind_set, molecule_set, ff_type, iw)
Determine the number of unique opbend kind and allocate opbend_kind_set based on the present improper...
subroutine, public force_field_pack_impr(particle_set, molecule_kind_set, molecule_set, ainfo, chm_info, inp_info, gro_info, iw)
Pack in impropers information needed for the force_field.
subroutine, public force_field_pack_radius(atomic_kind_set, iw, subsys_section)
Set up the radius of the electrostatic multipole in Fist.
subroutine, public force_field_pack_shell(particle_set, atomic_kind_set, molecule_kind_set, molecule_set, root_section, subsys_section, shell_particle_set, core_particle_set, cell, iw, inp_info)
Set up shell potential parameters.
subroutine, public force_field_pack_opbend(particle_set, molecule_kind_set, molecule_set, ainfo, inp_info, iw)
Pack in opbend information needed for the force_field. No loop over params for charmm,...
subroutine, public force_field_pack_charges(charges, charges_section, particle_set, my_qmmm, qmmm_env, inp_info, iw4)
Set up array of full charges.
subroutine, public force_field_pack_pol(atomic_kind_set, iw, inp_info)
Set up the polarizable FF parameters.
subroutine, public force_field_unique_bend(particle_set, molecule_kind_set, molecule_set, ff_type, iw)
Determine the number of unique bend kind and allocate bend_kind_set.
subroutine, public force_field_unique_impr(particle_set, molecule_kind_set, molecule_set, ff_type, iw)
Determine the number of unique impr kind and allocate impr_kind_set.
subroutine, public force_field_pack_splines(atomic_kind_set, ff_type, iw2, iw3, iw4, potparm, do_zbl, nonbonded_type)
create the pair potential spline environment
subroutine, public force_field_unique_tors(particle_set, molecule_kind_set, molecule_set, ff_type, iw)
Determine the number of unique torsion kind and allocate torsion_kind_set.
subroutine, public force_field_pack_nonbond(atomic_kind_set, ff_type, qmmm_env, fatal, iw, ainfo, chm_info, inp_info, gro_info, amb_info, potparm_nonbond, ewald_env)
Assign input and potential info to potparm_nonbond.
subroutine, public force_field_pack_charge(atomic_kind_set, qmmm_env, fatal, iw, iw4, ainfo, my_qmmm, inp_info)
Set up atomic_kind_set()fist_potential%[qeff] and shell potential parameters.
subroutine, public force_field_pack_bend(particle_set, molecule_kind_set, molecule_set, fatal, ainfo, chm_info, inp_info, gro_info, amb_info, iw)
Pack in bends information needed for the force_field.
subroutine, public force_field_pack_eicut(atomic_kind_set, ff_type, potparm_nonbond, ewald_env, iw)
Compute the electrostatic interaction cutoffs.
subroutine, public force_field_pack_bond(particle_set, molecule_kind_set, molecule_set, fatal, ainfo, chm_info, inp_info, gro_info, amb_info, iw)
Pack in bonds information needed for the force_field.
subroutine, public force_field_pack_damp(atomic_kind_set, iw, inp_info)
Set up damping parameters.
subroutine, public force_field_pack_nonbond14(atomic_kind_set, ff_type, qmmm_env, iw, ainfo, chm_info, inp_info, gro_info, amb_info, potparm_nonbond14, ewald_env)
Assign input and potential info to potparm_nonbond14.
subroutine, public force_field_pack_ub(particle_set, molecule_kind_set, molecule_set, ainfo, chm_info, inp_info, iw)
Pack in Urey-Bradley information needed for the force_field.
subroutine, public force_field_pack(particle_set, atomic_kind_set, molecule_kind_set, molecule_set, ewald_env, fist_nonbond_env, ff_type, root_section, qmmm, qmmm_env, mm_section, subsys_section, shell_particle_set, core_particle_set, cell)
Pack in all the information needed for the force_field.
subroutine, public get_generic_info(gen_section, func_name, xfunction, parameters, values, var_values, size_variables, i_rep_sec, input_variables)
Reads from the input structure all information for generic functions.
subroutine, public force_field_qeff_output(particle_set, molecule_kind_set, molecule_set, mm_section, charges)
Compute the total qeff charges for each molecule kind and total system.
subroutine, public clean_intra_force_kind(molecule_kind_set, mm_section)
Removes UNSET force field types.
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
integer, parameter, public default_path_length
Utility routines for the memory handling.
Define the molecule kind structure types and the corresponding functionality.
subroutine, public get_molecule_kind(molecule_kind, atom_list, bond_list, bend_list, ub_list, impr_list, opbend_list, colv_list, fixd_list, g3x3_list, g4x6_list, vsite_list, torsion_list, shell_list, name, mass, charge, kind_number, natom, nbend, nbond, nub, nimpr, nopbend, nconstraint, nconstraint_fixd, nfixd, ncolv, ng3x3, ng4x6, nvsite, nfixd_restraint, ng3x3_restraint, ng4x6_restraint, nvsite_restraint, nrestraints, nmolecule, nsgf, nshell, ntorsion, molecule_list, nelectron, nelectron_alpha, nelectron_beta, bond_kind_set, bend_kind_set, ub_kind_set, impr_kind_set, opbend_kind_set, torsion_kind_set, molname_generated)
Get informations about a molecule kind.
subroutine, public set_molecule_kind(molecule_kind, name, mass, charge, kind_number, molecule_list, atom_list, nbond, bond_list, nbend, bend_list, nub, ub_list, nimpr, impr_list, nopbend, opbend_list, ntorsion, torsion_list, fixd_list, ncolv, colv_list, ng3x3, g3x3_list, ng4x6, nfixd, g4x6_list, nvsite, vsite_list, ng3x3_restraint, ng4x6_restraint, nfixd_restraint, nshell, shell_list, nvsite_restraint, bond_kind_set, bend_kind_set, ub_kind_set, torsion_kind_set, impr_kind_set, opbend_kind_set, nelectron, nsgf, molname_generated)
Set the components of a molecule kind.
Define the data structure for the molecule information.
subroutine, public get_molecule(molecule, molecule_kind, lmi, lci, lg3x3, lg4x6, lcolv, first_atom, last_atom, first_shell, last_shell)
Get components from a molecule data set.
Define the data structure for the particle information.
Utilities for string manipulations.
subroutine, public compress(string, full)
Eliminate multiple space characters in a string. If full is .TRUE., then all spaces are eliminated.
Provides all information about an atomic kind.
Type defining parameters related to the simulation cell.
type of a logger, at the moment it contains just a print level starting at which level it should be l...
to build arrays of pointers