LCOV - code coverage report
Current view: top level - src/subsys - molecule_kind_types.F (source / functions) Hit Total Coverage
Test: CP2K Regtests (git:262480d) Lines: 291 302 96.4 %
Date: 2024-11-22 07:00:40 Functions: 8 24 33.3 %

          Line data    Source code
       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             : !> \brief Define the molecule kind structure types and the corresponding
      10             : !>      functionality
      11             : !> \par History
      12             : !>      Teodoro Laino [tlaino] 12.2008 - Preparing for VIRTUAL SITE constraints
      13             : !>                                       (patch by Marcel Baer)
      14             : !> \author MK (22.08.2003)
      15             : ! **************************************************************************************************
      16             : MODULE molecule_kind_types
      17             :    USE atomic_kind_types,               ONLY: atomic_kind_type,&
      18             :                                               get_atomic_kind
      19             :    USE cell_types,                      ONLY: use_perd_x,&
      20             :                                               use_perd_xy,&
      21             :                                               use_perd_xyz,&
      22             :                                               use_perd_xz,&
      23             :                                               use_perd_y,&
      24             :                                               use_perd_yz,&
      25             :                                               use_perd_z
      26             :    USE colvar_types,                    ONLY: &
      27             :         acid_hyd_dist_colvar_id, acid_hyd_shell_colvar_id, angle_colvar_id, colvar_counters, &
      28             :         combine_colvar_id, coord_colvar_id, dfunct_colvar_id, dist_colvar_id, gyration_colvar_id, &
      29             :         hydronium_dist_colvar_id, hydronium_shell_colvar_id, no_colvar_id, &
      30             :         plane_distance_colvar_id, plane_plane_angle_colvar_id, population_colvar_id, &
      31             :         qparm_colvar_id, reaction_path_colvar_id, rotation_colvar_id, torsion_colvar_id, &
      32             :         xyz_diag_colvar_id, xyz_outerdiag_colvar_id
      33             :    USE cp_log_handling,                 ONLY: cp_get_default_logger,&
      34             :                                               cp_logger_type
      35             :    USE cp_output_handling,              ONLY: cp_print_key_finished_output,&
      36             :                                               cp_print_key_unit_nr
      37             :    USE force_field_kind_types,          ONLY: &
      38             :         bend_kind_type, bond_kind_type, do_ff_undef, impr_kind_dealloc_ref, impr_kind_type, &
      39             :         opbend_kind_type, torsion_kind_dealloc_ref, torsion_kind_type, ub_kind_dealloc_ref, &
      40             :         ub_kind_type
      41             :    USE input_section_types,             ONLY: section_vals_type
      42             :    USE kinds,                           ONLY: default_string_length,&
      43             :                                               dp
      44             :    USE shell_potential_types,           ONLY: shell_kind_type
      45             : #include "../base/base_uses.f90"
      46             : 
      47             :    IMPLICIT NONE
      48             : 
      49             :    PRIVATE
      50             : 
      51             :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'molecule_kind_types'
      52             : 
      53             : ! *** Define the derived structure types ***
      54             : 
      55             : ! **************************************************************************************************
      56             :    TYPE atom_type
      57             :       TYPE(atomic_kind_type), POINTER :: atomic_kind => NULL()
      58             :       INTEGER :: id_name = 0
      59             :    END TYPE atom_type
      60             : 
      61             : ! **************************************************************************************************
      62             :    TYPE shell_type
      63             :       INTEGER :: a = 0
      64             :       CHARACTER(LEN=default_string_length)  :: name = ""
      65             :       TYPE(shell_kind_type), POINTER :: shell_kind => NULL()
      66             :    END TYPE shell_type
      67             : 
      68             : ! **************************************************************************************************
      69             :    TYPE bond_type
      70             :       INTEGER :: a = 0, b = 0
      71             :       INTEGER :: id_type = do_ff_undef, itype = 0
      72             :       TYPE(bond_kind_type), POINTER :: bond_kind => NULL()
      73             :    END TYPE bond_type
      74             : 
      75             : ! **************************************************************************************************
      76             :    TYPE bend_type
      77             :       INTEGER :: a = 0, b = 0, c = 0
      78             :       INTEGER :: id_type = do_ff_undef, itype = 0
      79             :       TYPE(bend_kind_type), POINTER :: bend_kind => NULL()
      80             :    END TYPE bend_type
      81             : 
      82             : ! **************************************************************************************************
      83             :    TYPE ub_type
      84             :       INTEGER :: a = 0, b = 0, c = 0
      85             :       INTEGER :: id_type = do_ff_undef, itype = 0
      86             :       TYPE(ub_kind_type), POINTER :: ub_kind => NULL()
      87             :    END TYPE ub_type
      88             : 
      89             : ! **************************************************************************************************
      90             :    TYPE torsion_type
      91             :       INTEGER :: a = 0, b = 0, c = 0, d = 0
      92             :       INTEGER :: id_type = do_ff_undef, itype = 0
      93             :       TYPE(torsion_kind_type), POINTER :: torsion_kind => NULL()
      94             :    END TYPE torsion_type
      95             : 
      96             : ! **************************************************************************************************
      97             :    TYPE impr_type
      98             :       INTEGER :: a = 0, b = 0, c = 0, d = 0
      99             :       INTEGER :: id_type = do_ff_undef, itype = 0
     100             :       TYPE(impr_kind_type), POINTER :: impr_kind => NULL()
     101             :    END TYPE impr_type
     102             : 
     103             : ! **************************************************************************************************
     104             :    TYPE opbend_type
     105             :       INTEGER :: a = 0, b = 0, c = 0, d = 0
     106             :       INTEGER :: id_type = do_ff_undef, itype = 0
     107             :       TYPE(opbend_kind_type), POINTER :: opbend_kind => NULL()
     108             :    END TYPE opbend_type
     109             : 
     110             : ! **************************************************************************************************
     111             :    TYPE restraint_type
     112             :       LOGICAL       :: active = .FALSE.
     113             :       REAL(KIND=dp) :: k0 = 0.0_dp
     114             :    END TYPE restraint_type
     115             : 
     116             : ! **************************************************************************************************
     117             :    TYPE colvar_constraint_type
     118             :       INTEGER                        :: type_id = no_colvar_id
     119             :       INTEGER                        :: inp_seq_num = 0
     120             :       LOGICAL                        :: use_points = .FALSE.
     121             :       REAL(KIND=dp)                :: expected_value = 0.0_dp
     122             :       REAL(KIND=dp)                :: expected_value_growth_speed = 0.0_dp
     123             :       INTEGER, POINTER, DIMENSION(:) :: i_atoms => NULL()
     124             :       TYPE(restraint_type)           :: restraint = restraint_type()
     125             :    END TYPE colvar_constraint_type
     126             : 
     127             : ! **************************************************************************************************
     128             :    TYPE g3x3_constraint_type
     129             :       INTEGER                        :: a = 0, b = 0, c = 0
     130             :       REAL(KIND=dp)                :: dab = 0.0_dp, dac = 0.0_dp, dbc = 0.0_dp
     131             :       TYPE(restraint_type)           :: restraint = restraint_type()
     132             :    END TYPE g3x3_constraint_type
     133             : 
     134             : ! **************************************************************************************************
     135             :    TYPE g4x6_constraint_type
     136             :       INTEGER                        :: a = 0, b = 0, c = 0, d = 0
     137             :       REAL(KIND=dp)                :: dab = 0.0_dp, dac = 0.0_dp, dbc = 0.0_dp, dad = 0.0_dp, dbd = 0.0_dp, dcd = 0.0_dp
     138             :       TYPE(restraint_type)           :: restraint = restraint_type()
     139             :    END TYPE g4x6_constraint_type
     140             : 
     141             : ! **************************************************************************************************
     142             :    TYPE vsite_constraint_type
     143             :       INTEGER                        :: a = 0, b = 0, c = 0, d = 0
     144             :       REAL(KIND=dp)                :: wbc = 0.0_dp, wdc = 0.0_dp
     145             :       TYPE(restraint_type)           :: restraint = restraint_type()
     146             :    END TYPE vsite_constraint_type
     147             : 
     148             : ! **************************************************************************************************
     149             :    TYPE fixd_constraint_type
     150             :       TYPE(restraint_type)           :: restraint = restraint_type()
     151             :       INTEGER                        :: fixd = 0, itype = 0
     152             :       REAL(KIND=dp), DIMENSION(3)    :: coord = 0.0_dp
     153             :    END TYPE fixd_constraint_type
     154             : 
     155             : ! **************************************************************************************************
     156             :    TYPE local_fixd_constraint_type
     157             :       INTEGER                        :: ifixd_index = 0, ikind = 0
     158             :    END TYPE local_fixd_constraint_type
     159             : 
     160             : ! **************************************************************************************************
     161             :    TYPE molecule_kind_type
     162             :       TYPE(atom_type), DIMENSION(:), POINTER            :: atom_list => NULL()
     163             :       TYPE(bond_kind_type), DIMENSION(:), POINTER       :: bond_kind_set => NULL()
     164             :       TYPE(bond_type), DIMENSION(:), POINTER            :: bond_list => NULL()
     165             :       TYPE(bend_kind_type), DIMENSION(:), POINTER       :: bend_kind_set => NULL()
     166             :       TYPE(bend_type), DIMENSION(:), POINTER            :: bend_list => NULL()
     167             :       TYPE(ub_kind_type), DIMENSION(:), POINTER         :: ub_kind_set => NULL()
     168             :       TYPE(ub_type), DIMENSION(:), POINTER              :: ub_list => NULL()
     169             :       TYPE(torsion_kind_type), DIMENSION(:), POINTER    :: torsion_kind_set => NULL()
     170             :       TYPE(torsion_type), DIMENSION(:), POINTER         :: torsion_list => NULL()
     171             :       TYPE(impr_kind_type), DIMENSION(:), POINTER       :: impr_kind_set => NULL()
     172             :       TYPE(impr_type), DIMENSION(:), POINTER            :: impr_list => NULL()
     173             :       TYPE(opbend_kind_type), DIMENSION(:), POINTER     :: opbend_kind_set => NULL()
     174             :       TYPE(opbend_type), DIMENSION(:), POINTER          :: opbend_list => NULL()
     175             :       TYPE(colvar_constraint_type), DIMENSION(:), POINTER :: colv_list => NULL()
     176             :       TYPE(g3x3_constraint_type), DIMENSION(:), POINTER   :: g3x3_list => NULL()
     177             :       TYPE(g4x6_constraint_type), DIMENSION(:), POINTER   :: g4x6_list => NULL()
     178             :       TYPE(vsite_constraint_type), DIMENSION(:), POINTER  :: vsite_list => NULL()
     179             :       TYPE(fixd_constraint_type), DIMENSION(:), POINTER   :: fixd_list => NULL()
     180             :       TYPE(shell_type), DIMENSION(:), POINTER           :: shell_list => NULL()
     181             :       CHARACTER(LEN=default_string_length)              :: name = ""
     182             :       REAL(KIND=dp)                                   :: charge = 0.0_dp, &
     183             :                                                          mass = 0.0_dp
     184             :       INTEGER                                           :: kind_number = 0, &
     185             :                                                            natom = 0, &
     186             :                                                            nbond = 0, &
     187             :                                                            nbend = 0, &
     188             :                                                            nimpr = 0, &
     189             :                                                            nopbend = 0, &
     190             :                                                            ntorsion = 0, &
     191             :                                                            nub = 0, &
     192             :                                                            ng3x3 = 0, ng3x3_restraint = 0, &
     193             :                                                            ng4x6 = 0, ng4x6_restraint = 0, &
     194             :                                                            nvsite = 0, nvsite_restraint = 0, &
     195             :                                                            nfixd = 0, nfixd_restraint = 0, &
     196             :                                                            nmolecule = 0, nshell = 0
     197             :       TYPE(colvar_counters)                             :: ncolv = colvar_counters()
     198             :       INTEGER                                           :: nsgf = 0, nelectron = 0, &
     199             :                                                            nelectron_alpha = 0, &
     200             :                                                            nelectron_beta = 0
     201             :       INTEGER, DIMENSION(:), POINTER                    :: molecule_list => NULL()
     202             :       LOGICAL                                           :: molname_generated = .FALSE.
     203             :    END TYPE molecule_kind_type
     204             : 
     205             :    ! *** Public subroutines ***
     206             :    PUBLIC :: allocate_molecule_kind_set, &
     207             :              deallocate_molecule_kind_set, &
     208             :              get_molecule_kind, &
     209             :              get_molecule_kind_set, &
     210             :              set_molecule_kind, &
     211             :              write_molecule_kind_set, &
     212             :              setup_colvar_counters
     213             : 
     214             :    ! *** Public data types ***
     215             :    PUBLIC :: atom_type, &
     216             :              bend_type, &
     217             :              bond_type, &
     218             :              ub_type, &
     219             :              torsion_type, &
     220             :              impr_type, &
     221             :              opbend_type, &
     222             :              colvar_constraint_type, &
     223             :              g3x3_constraint_type, &
     224             :              g4x6_constraint_type, &
     225             :              vsite_constraint_type, &
     226             :              fixd_constraint_type, &
     227             :              local_fixd_constraint_type, &
     228             :              molecule_kind_type, &
     229             :              shell_type
     230             : 
     231             : CONTAINS
     232             : 
     233             : ! **************************************************************************************************
     234             : !> \brief ...
     235             : !> \param colv_list ...
     236             : !> \param ncolv ...
     237             : ! **************************************************************************************************
     238      142027 :    SUBROUTINE setup_colvar_counters(colv_list, ncolv)
     239             :       TYPE(colvar_constraint_type), DIMENSION(:), &
     240             :          POINTER                                         :: colv_list
     241             :       TYPE(colvar_counters), INTENT(OUT)                 :: ncolv
     242             : 
     243             :       INTEGER                                            :: k
     244             : 
     245      142027 :       IF (ASSOCIATED(colv_list)) THEN
     246        1070 :          DO k = 1, SIZE(colv_list)
     247         448 :             IF (colv_list(k)%restraint%active) ncolv%nrestraint = ncolv%nrestraint + 1
     248         622 :             SELECT CASE (colv_list(k)%type_id)
     249             :             CASE (angle_colvar_id)
     250          50 :                ncolv%nangle = ncolv%nangle + 1
     251             :             CASE (coord_colvar_id)
     252           2 :                ncolv%ncoord = ncolv%ncoord + 1
     253             :             CASE (population_colvar_id)
     254           0 :                ncolv%npopulation = ncolv%npopulation + 1
     255             :             CASE (gyration_colvar_id)
     256           0 :                ncolv%ngyration = ncolv%ngyration + 1
     257             :             CASE (rotation_colvar_id)
     258           0 :                ncolv%nrot = ncolv%nrot + 1
     259             :             CASE (dist_colvar_id)
     260         334 :                ncolv%ndist = ncolv%ndist + 1
     261             :             CASE (dfunct_colvar_id)
     262           4 :                ncolv%ndfunct = ncolv%ndfunct + 1
     263             :             CASE (plane_distance_colvar_id)
     264           0 :                ncolv%nplane_dist = ncolv%nplane_dist + 1
     265             :             CASE (plane_plane_angle_colvar_id)
     266           4 :                ncolv%nplane_angle = ncolv%nplane_angle + 1
     267             :             CASE (torsion_colvar_id)
     268          38 :                ncolv%ntorsion = ncolv%ntorsion + 1
     269             :             CASE (qparm_colvar_id)
     270           0 :                ncolv%nqparm = ncolv%nqparm + 1
     271             :             CASE (xyz_diag_colvar_id)
     272           6 :                ncolv%nxyz_diag = ncolv%nxyz_diag + 1
     273             :             CASE (xyz_outerdiag_colvar_id)
     274           6 :                ncolv%nxyz_outerdiag = ncolv%nxyz_outerdiag + 1
     275             :             CASE (hydronium_shell_colvar_id)
     276           0 :                ncolv%nhydronium_shell = ncolv%nhydronium_shell + 1
     277             :             CASE (hydronium_dist_colvar_id)
     278           0 :                ncolv%nhydronium_dist = ncolv%nhydronium_dist + 1
     279             :             CASE (acid_hyd_dist_colvar_id)
     280           0 :                ncolv%nacid_hyd_dist = ncolv%nacid_hyd_dist + 1
     281             :             CASE (acid_hyd_shell_colvar_id)
     282           0 :                ncolv%nacid_hyd_shell = ncolv%nacid_hyd_shell + 1
     283             :             CASE (reaction_path_colvar_id)
     284           2 :                ncolv%nreactionpath = ncolv%nreactionpath + 1
     285             :             CASE (combine_colvar_id)
     286           2 :                ncolv%ncombinecvs = ncolv%ncombinecvs + 1
     287             :             CASE DEFAULT
     288         448 :                CPABORT("")
     289             :             END SELECT
     290             :          END DO
     291             :       END IF
     292             :       ncolv%ntot = ncolv%ndist + &
     293             :                    ncolv%nangle + &
     294             :                    ncolv%ntorsion + &
     295             :                    ncolv%ncoord + &
     296             :                    ncolv%nplane_dist + &
     297             :                    ncolv%nplane_angle + &
     298             :                    ncolv%ndfunct + &
     299             :                    ncolv%nrot + &
     300             :                    ncolv%nqparm + &
     301             :                    ncolv%nxyz_diag + &
     302             :                    ncolv%nxyz_outerdiag + &
     303             :                    ncolv%nhydronium_shell + &
     304             :                    ncolv%nhydronium_dist + &
     305             :                    ncolv%nacid_hyd_dist + &
     306             :                    ncolv%nacid_hyd_shell + &
     307             :                    ncolv%nreactionpath + &
     308             :                    ncolv%ncombinecvs + &
     309             :                    ncolv%npopulation + &
     310      142027 :                    ncolv%ngyration
     311             : 
     312      142027 :    END SUBROUTINE setup_colvar_counters
     313             : 
     314             : ! **************************************************************************************************
     315             : !> \brief   Allocate and initialize a molecule kind set.
     316             : !> \param molecule_kind_set ...
     317             : !> \param nmolecule_kind ...
     318             : !> \date    22.08.2003
     319             : !> \author  MK
     320             : !> \version 1.0
     321             : ! **************************************************************************************************
     322        9512 :    SUBROUTINE allocate_molecule_kind_set(molecule_kind_set, nmolecule_kind)
     323             :       TYPE(molecule_kind_type), DIMENSION(:), POINTER    :: molecule_kind_set
     324             :       INTEGER, INTENT(IN)                                :: nmolecule_kind
     325             : 
     326             :       INTEGER                                            :: imolecule_kind
     327             : 
     328        9512 :       IF (ASSOCIATED(molecule_kind_set)) THEN
     329           0 :          CALL deallocate_molecule_kind_set(molecule_kind_set)
     330             :       END IF
     331             : 
     332      160957 :       ALLOCATE (molecule_kind_set(nmolecule_kind))
     333             : 
     334      141933 :       DO imolecule_kind = 1, nmolecule_kind
     335      132421 :          molecule_kind_set(imolecule_kind)%kind_number = imolecule_kind
     336             :          CALL setup_colvar_counters(molecule_kind_set(imolecule_kind)%colv_list, &
     337      141933 :                                     molecule_kind_set(imolecule_kind)%ncolv)
     338             :       END DO
     339             : 
     340        9512 :    END SUBROUTINE allocate_molecule_kind_set
     341             : 
     342             : ! **************************************************************************************************
     343             : !> \brief   Deallocate a molecule kind set.
     344             : !> \param molecule_kind_set ...
     345             : !> \date    22.08.2003
     346             : !> \author  MK
     347             : !> \version 1.0
     348             : ! **************************************************************************************************
     349        9512 :    SUBROUTINE deallocate_molecule_kind_set(molecule_kind_set)
     350             : 
     351             :       TYPE(molecule_kind_type), DIMENSION(:), POINTER    :: molecule_kind_set
     352             : 
     353             :       INTEGER                                            :: i, imolecule_kind, j, nmolecule_kind
     354             : 
     355        9512 :       IF (ASSOCIATED(molecule_kind_set)) THEN
     356             : 
     357        9512 :          nmolecule_kind = SIZE(molecule_kind_set)
     358             : 
     359      141933 :          DO imolecule_kind = 1, nmolecule_kind
     360             : 
     361      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%atom_list)) THEN
     362      132421 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%atom_list)
     363             :             END IF
     364      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%bend_kind_set)) THEN
     365      122739 :                DO i = 1, SIZE(molecule_kind_set(imolecule_kind)%bend_kind_set)
     366       93640 :                   IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%bend_kind_set(i)%legendre%coeffs)) &
     367       31168 :                      DEALLOCATE (molecule_kind_set(imolecule_kind)%bend_kind_set(i)%legendre%coeffs)
     368             :                END DO
     369       29099 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%bend_kind_set)
     370             :             END IF
     371      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%bend_list)) THEN
     372      132421 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%bend_list)
     373             :             END IF
     374      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%ub_list)) THEN
     375      132421 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%ub_list)
     376             :             END IF
     377      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%ub_kind_set)) THEN
     378       29085 :                CALL ub_kind_dealloc_ref(molecule_kind_set(imolecule_kind)%ub_kind_set)
     379             :             END IF
     380      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%impr_list)) THEN
     381      132421 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%impr_list)
     382             :             END IF
     383      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%impr_kind_set)) THEN
     384        4882 :                DO i = 1, SIZE(molecule_kind_set(imolecule_kind)%impr_kind_set)
     385        4882 :                   CALL impr_kind_dealloc_ref() !This Subroutine doesn't deallocate anything, maybe needs to be implemented
     386             :                END DO
     387        1672 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%impr_kind_set)
     388             :             END IF
     389      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%opbend_list)) THEN
     390      132421 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%opbend_list)
     391             :             END IF
     392      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%opbend_kind_set)) THEN
     393        1672 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%opbend_kind_set)
     394             :             END IF
     395      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%bond_kind_set)) THEN
     396       29429 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%bond_kind_set)
     397             :             END IF
     398      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%bond_list)) THEN
     399      132421 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%bond_list)
     400             :             END IF
     401      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%colv_list)) THEN
     402         960 :                DO j = 1, SIZE(molecule_kind_set(imolecule_kind)%colv_list)
     403         960 :                   DEALLOCATE (molecule_kind_set(imolecule_kind)%colv_list(j)%i_atoms)
     404             :                END DO
     405         578 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%colv_list)
     406             :             END IF
     407      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%g3x3_list)) THEN
     408         270 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%g3x3_list)
     409             :             END IF
     410      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%g4x6_list)) THEN
     411          20 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%g4x6_list)
     412             :             END IF
     413      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%vsite_list)) THEN
     414          10 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%vsite_list)
     415             :             END IF
     416      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%fixd_list)) THEN
     417        4926 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%fixd_list)
     418             :             END IF
     419      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%torsion_kind_set)) THEN
     420       83961 :                DO i = 1, SIZE(molecule_kind_set(imolecule_kind)%torsion_kind_set)
     421       83961 :                   CALL torsion_kind_dealloc_ref(molecule_kind_set(imolecule_kind)%torsion_kind_set(i))
     422             :                END DO
     423        5532 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%torsion_kind_set)
     424             :             END IF
     425      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%shell_list)) THEN
     426       10876 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%shell_list)
     427             :             END IF
     428      132421 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%torsion_list)) THEN
     429      132421 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%torsion_list)
     430             :             END IF
     431      141933 :             IF (ASSOCIATED(molecule_kind_set(imolecule_kind)%molecule_list)) THEN
     432      132421 :                DEALLOCATE (molecule_kind_set(imolecule_kind)%molecule_list)
     433             :             END IF
     434             :          END DO
     435             : 
     436        9512 :          DEALLOCATE (molecule_kind_set)
     437             :       END IF
     438        9512 :       NULLIFY (molecule_kind_set)
     439             : 
     440        9512 :    END SUBROUTINE deallocate_molecule_kind_set
     441             : 
     442             : ! **************************************************************************************************
     443             : !> \brief   Get informations about a molecule kind.
     444             : !> \param molecule_kind ...
     445             : !> \param atom_list ...
     446             : !> \param bond_list ...
     447             : !> \param bend_list ...
     448             : !> \param ub_list ...
     449             : !> \param impr_list ...
     450             : !> \param opbend_list ...
     451             : !> \param colv_list ...
     452             : !> \param fixd_list ...
     453             : !> \param g3x3_list ...
     454             : !> \param g4x6_list ...
     455             : !> \param vsite_list ...
     456             : !> \param torsion_list ...
     457             : !> \param shell_list ...
     458             : !> \param name ...
     459             : !> \param mass ...
     460             : !> \param charge ...
     461             : !> \param kind_number ...
     462             : !> \param natom ...
     463             : !> \param nbend ...
     464             : !> \param nbond ...
     465             : !> \param nub ...
     466             : !> \param nimpr ...
     467             : !> \param nopbend ...
     468             : !> \param nconstraint ...
     469             : !> \param nconstraint_fixd ...
     470             : !> \param nfixd ...
     471             : !> \param ncolv ...
     472             : !> \param ng3x3 ...
     473             : !> \param ng4x6 ...
     474             : !> \param nvsite ...
     475             : !> \param nfixd_restraint ...
     476             : !> \param ng3x3_restraint ...
     477             : !> \param ng4x6_restraint ...
     478             : !> \param nvsite_restraint ...
     479             : !> \param nrestraints ...
     480             : !> \param nmolecule ...
     481             : !> \param nsgf ...
     482             : !> \param nshell ...
     483             : !> \param ntorsion ...
     484             : !> \param molecule_list ...
     485             : !> \param nelectron ...
     486             : !> \param nelectron_alpha ...
     487             : !> \param nelectron_beta ...
     488             : !> \param bond_kind_set ...
     489             : !> \param bend_kind_set ...
     490             : !> \param ub_kind_set ...
     491             : !> \param impr_kind_set ...
     492             : !> \param opbend_kind_set ...
     493             : !> \param torsion_kind_set ...
     494             : !> \param molname_generated ...
     495             : !> \date    27.08.2003
     496             : !> \author  MK
     497             : !> \version 1.0
     498             : ! **************************************************************************************************
     499    16345008 :    SUBROUTINE get_molecule_kind(molecule_kind, atom_list, bond_list, bend_list, &
     500             :                                 ub_list, impr_list, opbend_list, colv_list, fixd_list, &
     501             :                                 g3x3_list, g4x6_list, vsite_list, torsion_list, shell_list, &
     502             :                                 name, mass, charge, kind_number, natom, nbend, nbond, nub, &
     503             :                                 nimpr, nopbend, nconstraint, nconstraint_fixd, nfixd, ncolv, ng3x3, ng4x6, &
     504             :                                 nvsite, nfixd_restraint, ng3x3_restraint, ng4x6_restraint, &
     505             :                                 nvsite_restraint, nrestraints, nmolecule, nsgf, nshell, ntorsion, &
     506             :                                 molecule_list, nelectron, nelectron_alpha, nelectron_beta, &
     507             :                                 bond_kind_set, bend_kind_set, &
     508             :                                 ub_kind_set, impr_kind_set, opbend_kind_set, torsion_kind_set, &
     509             :                                 molname_generated)
     510             : 
     511             :       TYPE(molecule_kind_type), INTENT(IN)               :: molecule_kind
     512             :       TYPE(atom_type), DIMENSION(:), OPTIONAL, POINTER   :: atom_list
     513             :       TYPE(bond_type), DIMENSION(:), OPTIONAL, POINTER   :: bond_list
     514             :       TYPE(bend_type), DIMENSION(:), OPTIONAL, POINTER   :: bend_list
     515             :       TYPE(ub_type), DIMENSION(:), OPTIONAL, POINTER     :: ub_list
     516             :       TYPE(impr_type), DIMENSION(:), OPTIONAL, POINTER   :: impr_list
     517             :       TYPE(opbend_type), DIMENSION(:), OPTIONAL, POINTER :: opbend_list
     518             :       TYPE(colvar_constraint_type), DIMENSION(:), &
     519             :          OPTIONAL, POINTER                               :: colv_list
     520             :       TYPE(fixd_constraint_type), DIMENSION(:), &
     521             :          OPTIONAL, POINTER                               :: fixd_list
     522             :       TYPE(g3x3_constraint_type), DIMENSION(:), &
     523             :          OPTIONAL, POINTER                               :: g3x3_list
     524             :       TYPE(g4x6_constraint_type), DIMENSION(:), &
     525             :          OPTIONAL, POINTER                               :: g4x6_list
     526             :       TYPE(vsite_constraint_type), DIMENSION(:), &
     527             :          OPTIONAL, POINTER                               :: vsite_list
     528             :       TYPE(torsion_type), DIMENSION(:), OPTIONAL, &
     529             :          POINTER                                         :: torsion_list
     530             :       TYPE(shell_type), DIMENSION(:), OPTIONAL, POINTER  :: shell_list
     531             :       CHARACTER(LEN=default_string_length), &
     532             :          INTENT(OUT), OPTIONAL                           :: name
     533             :       REAL(KIND=dp), OPTIONAL                            :: mass, charge
     534             :       INTEGER, INTENT(OUT), OPTIONAL                     :: kind_number, natom, nbend, nbond, nub, &
     535             :                                                             nimpr, nopbend, nconstraint, &
     536             :                                                             nconstraint_fixd, nfixd
     537             :       TYPE(colvar_counters), INTENT(out), OPTIONAL       :: ncolv
     538             :       INTEGER, INTENT(OUT), OPTIONAL :: ng3x3, ng4x6, nvsite, nfixd_restraint, ng3x3_restraint, &
     539             :          ng4x6_restraint, nvsite_restraint, nrestraints, nmolecule, nsgf, nshell, ntorsion
     540             :       INTEGER, DIMENSION(:), OPTIONAL, POINTER           :: molecule_list
     541             :       INTEGER, INTENT(OUT), OPTIONAL                     :: nelectron, nelectron_alpha, &
     542             :                                                             nelectron_beta
     543             :       TYPE(bond_kind_type), DIMENSION(:), OPTIONAL, &
     544             :          POINTER                                         :: bond_kind_set
     545             :       TYPE(bend_kind_type), DIMENSION(:), OPTIONAL, &
     546             :          POINTER                                         :: bend_kind_set
     547             :       TYPE(ub_kind_type), DIMENSION(:), OPTIONAL, &
     548             :          POINTER                                         :: ub_kind_set
     549             :       TYPE(impr_kind_type), DIMENSION(:), OPTIONAL, &
     550             :          POINTER                                         :: impr_kind_set
     551             :       TYPE(opbend_kind_type), DIMENSION(:), OPTIONAL, &
     552             :          POINTER                                         :: opbend_kind_set
     553             :       TYPE(torsion_kind_type), DIMENSION(:), OPTIONAL, &
     554             :          POINTER                                         :: torsion_kind_set
     555             :       LOGICAL, INTENT(OUT), OPTIONAL                     :: molname_generated
     556             : 
     557             :       INTEGER                                            :: i
     558             : 
     559    16345008 :       IF (PRESENT(atom_list)) atom_list => molecule_kind%atom_list
     560    16345008 :       IF (PRESENT(bend_list)) bend_list => molecule_kind%bend_list
     561    16345008 :       IF (PRESENT(bond_list)) bond_list => molecule_kind%bond_list
     562    16345008 :       IF (PRESENT(impr_list)) impr_list => molecule_kind%impr_list
     563    16345008 :       IF (PRESENT(opbend_list)) opbend_list => molecule_kind%opbend_list
     564    16345008 :       IF (PRESENT(ub_list)) ub_list => molecule_kind%ub_list
     565    16345008 :       IF (PRESENT(bond_kind_set)) bond_kind_set => molecule_kind%bond_kind_set
     566    16345008 :       IF (PRESENT(bend_kind_set)) bend_kind_set => molecule_kind%bend_kind_set
     567    16345008 :       IF (PRESENT(ub_kind_set)) ub_kind_set => molecule_kind%ub_kind_set
     568    16345008 :       IF (PRESENT(impr_kind_set)) impr_kind_set => molecule_kind%impr_kind_set
     569    16345008 :       IF (PRESENT(opbend_kind_set)) opbend_kind_set => molecule_kind%opbend_kind_set
     570    16345008 :       IF (PRESENT(torsion_kind_set)) torsion_kind_set => molecule_kind%torsion_kind_set
     571    16345008 :       IF (PRESENT(colv_list)) colv_list => molecule_kind%colv_list
     572    16345008 :       IF (PRESENT(g3x3_list)) g3x3_list => molecule_kind%g3x3_list
     573    16345008 :       IF (PRESENT(g4x6_list)) g4x6_list => molecule_kind%g4x6_list
     574    16345008 :       IF (PRESENT(vsite_list)) vsite_list => molecule_kind%vsite_list
     575    16345008 :       IF (PRESENT(fixd_list)) fixd_list => molecule_kind%fixd_list
     576    16345008 :       IF (PRESENT(torsion_list)) torsion_list => molecule_kind%torsion_list
     577    16345008 :       IF (PRESENT(shell_list)) shell_list => molecule_kind%shell_list
     578    16345008 :       IF (PRESENT(name)) name = molecule_kind%name
     579    16345008 :       IF (PRESENT(molname_generated)) molname_generated = molecule_kind%molname_generated
     580    16345008 :       IF (PRESENT(mass)) mass = molecule_kind%mass
     581    16345008 :       IF (PRESENT(charge)) charge = molecule_kind%charge
     582    16345008 :       IF (PRESENT(kind_number)) kind_number = molecule_kind%kind_number
     583    16345008 :       IF (PRESENT(natom)) natom = molecule_kind%natom
     584    16345008 :       IF (PRESENT(nbend)) nbend = molecule_kind%nbend
     585    16345008 :       IF (PRESENT(nbond)) nbond = molecule_kind%nbond
     586    16345008 :       IF (PRESENT(nub)) nub = molecule_kind%nub
     587    16345008 :       IF (PRESENT(nimpr)) nimpr = molecule_kind%nimpr
     588    16345008 :       IF (PRESENT(nopbend)) nopbend = molecule_kind%nopbend
     589    16345008 :       IF (PRESENT(nconstraint)) nconstraint = (molecule_kind%ncolv%ntot - molecule_kind%ncolv%nrestraint) + &
     590             :                                               3*(molecule_kind%ng3x3 - molecule_kind%ng3x3_restraint) + &
     591             :                                               6*(molecule_kind%ng4x6 - molecule_kind%ng4x6_restraint) + &
     592     3172209 :                                               3*(molecule_kind%nvsite - molecule_kind%nvsite_restraint)
     593    16345008 :       IF (PRESENT(ncolv)) ncolv = molecule_kind%ncolv
     594    16345008 :       IF (PRESENT(ng3x3)) ng3x3 = molecule_kind%ng3x3
     595    16345008 :       IF (PRESENT(ng4x6)) ng4x6 = molecule_kind%ng4x6
     596    16345008 :       IF (PRESENT(nvsite)) nvsite = molecule_kind%nvsite
     597             :       ! Number of atoms that have one or more components fixed
     598    16345008 :       IF (PRESENT(nfixd)) nfixd = molecule_kind%nfixd
     599             :       ! Number of degrees of freedom fixed
     600    16345008 :       IF (PRESENT(nconstraint_fixd)) THEN
     601      275670 :          nconstraint_fixd = 0
     602      275670 :          IF (molecule_kind%nfixd /= 0) THEN
     603      171910 :             DO i = 1, SIZE(molecule_kind%fixd_list)
     604      170016 :                IF (molecule_kind%fixd_list(i)%restraint%active) CYCLE
     605        1894 :                SELECT CASE (molecule_kind%fixd_list(i)%itype)
     606             :                CASE (use_perd_x, use_perd_y, use_perd_z)
     607       62976 :                   nconstraint_fixd = nconstraint_fixd + 1
     608             :                CASE (use_perd_xy, use_perd_xz, use_perd_yz)
     609       20992 :                   nconstraint_fixd = nconstraint_fixd + 2
     610             :                CASE (use_perd_xyz)
     611      169588 :                   nconstraint_fixd = nconstraint_fixd + 3
     612             :                END SELECT
     613             :             END DO
     614             :          END IF
     615             :       END IF
     616    16345008 :       IF (PRESENT(ng3x3_restraint)) ng3x3_restraint = molecule_kind%ng3x3_restraint
     617    16345008 :       IF (PRESENT(ng4x6_restraint)) ng4x6_restraint = molecule_kind%ng4x6_restraint
     618    16345008 :       IF (PRESENT(nvsite_restraint)) nvsite_restraint = molecule_kind%nvsite_restraint
     619    16345008 :       IF (PRESENT(nfixd_restraint)) nfixd_restraint = molecule_kind%nfixd_restraint
     620    16345008 :       IF (PRESENT(nrestraints)) nrestraints = molecule_kind%ncolv%nrestraint + &
     621             :                                               molecule_kind%ng3x3_restraint + &
     622             :                                               molecule_kind%ng4x6_restraint + &
     623      261086 :                                               molecule_kind%nvsite_restraint
     624    16345008 :       IF (PRESENT(nmolecule)) nmolecule = molecule_kind%nmolecule
     625    16345008 :       IF (PRESENT(nshell)) nshell = molecule_kind%nshell
     626    16345008 :       IF (PRESENT(ntorsion)) ntorsion = molecule_kind%ntorsion
     627    16345008 :       IF (PRESENT(nsgf)) nsgf = molecule_kind%nsgf
     628    16345008 :       IF (PRESENT(nelectron)) nelectron = molecule_kind%nelectron
     629    16345008 :       IF (PRESENT(nelectron_alpha)) nelectron_alpha = molecule_kind%nelectron_beta
     630    16345008 :       IF (PRESENT(nelectron_beta)) nelectron_beta = molecule_kind%nelectron_alpha
     631    16345008 :       IF (PRESENT(molecule_list)) molecule_list => molecule_kind%molecule_list
     632             : 
     633    16345008 :    END SUBROUTINE get_molecule_kind
     634             : 
     635             : ! **************************************************************************************************
     636             : !> \brief   Get informations about a molecule kind set.
     637             : !> \param molecule_kind_set ...
     638             : !> \param maxatom ...
     639             : !> \param natom ...
     640             : !> \param nbond ...
     641             : !> \param nbend ...
     642             : !> \param nub ...
     643             : !> \param ntorsion ...
     644             : !> \param nimpr ...
     645             : !> \param nopbend ...
     646             : !> \param nconstraint ...
     647             : !> \param nconstraint_fixd ...
     648             : !> \param nmolecule ...
     649             : !> \param nrestraints ...
     650             : !> \date    27.08.2003
     651             : !> \author  MK
     652             : !> \version 1.0
     653             : ! **************************************************************************************************
     654       48228 :    SUBROUTINE get_molecule_kind_set(molecule_kind_set, maxatom, natom, &
     655             :                                     nbond, nbend, nub, ntorsion, nimpr, nopbend, &
     656             :                                     nconstraint, nconstraint_fixd, nmolecule, &
     657             :                                     nrestraints)
     658             : 
     659             :       TYPE(molecule_kind_type), DIMENSION(:), INTENT(IN) :: molecule_kind_set
     660             :       INTEGER, INTENT(OUT), OPTIONAL                     :: maxatom, natom, nbond, nbend, nub, &
     661             :                                                             ntorsion, nimpr, nopbend, nconstraint, &
     662             :                                                             nconstraint_fixd, nmolecule, &
     663             :                                                             nrestraints
     664             : 
     665             :       INTEGER :: ibend, ibond, iimpr, imolecule_kind, iopbend, itorsion, iub, na, nc, nc_fixd, &
     666             :          nfixd_restraint, nm, nmolecule_kind, nrestraints_tot
     667             : 
     668       48228 :       IF (PRESENT(maxatom)) maxatom = 0
     669       48228 :       IF (PRESENT(natom)) natom = 0
     670       48228 :       IF (PRESENT(nbond)) nbond = 0
     671       48228 :       IF (PRESENT(nbend)) nbend = 0
     672       48228 :       IF (PRESENT(nub)) nub = 0
     673       48228 :       IF (PRESENT(ntorsion)) ntorsion = 0
     674       48228 :       IF (PRESENT(nimpr)) nimpr = 0
     675       48228 :       IF (PRESENT(nopbend)) nopbend = 0
     676       48228 :       IF (PRESENT(nconstraint)) nconstraint = 0
     677       48228 :       IF (PRESENT(nconstraint_fixd)) nconstraint_fixd = 0
     678       48228 :       IF (PRESENT(nrestraints)) nrestraints = 0
     679       48228 :       IF (PRESENT(nmolecule)) nmolecule = 0
     680             : 
     681       48228 :       nmolecule_kind = SIZE(molecule_kind_set)
     682             : 
     683      309314 :       DO imolecule_kind = 1, nmolecule_kind
     684       48228 :          ASSOCIATE (molecule_kind => molecule_kind_set(imolecule_kind))
     685             : 
     686             :             CALL get_molecule_kind(molecule_kind=molecule_kind, &
     687             :                                    natom=na, &
     688             :                                    nbond=ibond, &
     689             :                                    nbend=ibend, &
     690             :                                    nub=iub, &
     691             :                                    ntorsion=itorsion, &
     692             :                                    nimpr=iimpr, &
     693             :                                    nopbend=iopbend, &
     694             :                                    nconstraint=nc, &
     695             :                                    nconstraint_fixd=nc_fixd, &
     696             :                                    nfixd_restraint=nfixd_restraint, &
     697             :                                    nrestraints=nrestraints_tot, &
     698      261086 :                                    nmolecule=nm)
     699      261086 :             IF (PRESENT(maxatom)) maxatom = MAX(maxatom, na)
     700      261086 :             IF (PRESENT(natom)) natom = natom + na*nm
     701      261086 :             IF (PRESENT(nbond)) nbond = nbond + ibond*nm
     702      261086 :             IF (PRESENT(nbend)) nbend = nbend + ibend*nm
     703      261086 :             IF (PRESENT(nub)) nub = nub + iub*nm
     704      261086 :             IF (PRESENT(ntorsion)) ntorsion = ntorsion + itorsion*nm
     705      261086 :             IF (PRESENT(nimpr)) nimpr = nimpr + iimpr*nm
     706      261086 :             IF (PRESENT(nopbend)) nopbend = nopbend + iopbend*nm
     707      261086 :             IF (PRESENT(nconstraint)) nconstraint = nconstraint + nc*nm + nc_fixd
     708      261086 :             IF (PRESENT(nconstraint_fixd)) nconstraint_fixd = nconstraint_fixd + nc_fixd
     709      261086 :             IF (PRESENT(nmolecule)) nmolecule = nmolecule + nm
     710      522172 :             IF (PRESENT(nrestraints)) nrestraints = nrestraints + nm*nrestraints_tot + nfixd_restraint
     711             : 
     712             :          END ASSOCIATE
     713             :       END DO
     714             : 
     715       48228 :    END SUBROUTINE get_molecule_kind_set
     716             : 
     717             : ! **************************************************************************************************
     718             : !> \brief   Set the components of a molecule kind.
     719             : !> \param molecule_kind ...
     720             : !> \param name ...
     721             : !> \param mass ...
     722             : !> \param charge ...
     723             : !> \param kind_number ...
     724             : !> \param molecule_list ...
     725             : !> \param atom_list ...
     726             : !> \param nbond ...
     727             : !> \param bond_list ...
     728             : !> \param nbend ...
     729             : !> \param bend_list ...
     730             : !> \param nub ...
     731             : !> \param ub_list ...
     732             : !> \param nimpr ...
     733             : !> \param impr_list ...
     734             : !> \param nopbend ...
     735             : !> \param opbend_list ...
     736             : !> \param ntorsion ...
     737             : !> \param torsion_list ...
     738             : !> \param fixd_list ...
     739             : !> \param ncolv ...
     740             : !> \param colv_list ...
     741             : !> \param ng3x3 ...
     742             : !> \param g3x3_list ...
     743             : !> \param ng4x6 ...
     744             : !> \param nfixd ...
     745             : !> \param g4x6_list ...
     746             : !> \param nvsite ...
     747             : !> \param vsite_list ...
     748             : !> \param ng3x3_restraint ...
     749             : !> \param ng4x6_restraint ...
     750             : !> \param nfixd_restraint ...
     751             : !> \param nshell ...
     752             : !> \param shell_list ...
     753             : !> \param nvsite_restraint ...
     754             : !> \param bond_kind_set ...
     755             : !> \param bend_kind_set ...
     756             : !> \param ub_kind_set ...
     757             : !> \param torsion_kind_set ...
     758             : !> \param impr_kind_set ...
     759             : !> \param opbend_kind_set ...
     760             : !> \param nelectron ...
     761             : !> \param nsgf ...
     762             : !> \param molname_generated ...
     763             : !> \date    27.08.2003
     764             : !> \author  MK
     765             : !> \version 1.0
     766             : ! **************************************************************************************************
     767     2005801 :    SUBROUTINE set_molecule_kind(molecule_kind, name, mass, charge, kind_number, &
     768             :                                 molecule_list, atom_list, nbond, bond_list, &
     769             :                                 nbend, bend_list, nub, ub_list, nimpr, impr_list, &
     770             :                                 nopbend, opbend_list, ntorsion, &
     771             :                                 torsion_list, fixd_list, ncolv, colv_list, ng3x3, &
     772             :                                 g3x3_list, ng4x6, nfixd, g4x6_list, nvsite, &
     773             :                                 vsite_list, ng3x3_restraint, ng4x6_restraint, &
     774             :                                 nfixd_restraint, nshell, shell_list, &
     775             :                                 nvsite_restraint, bond_kind_set, bend_kind_set, &
     776             :                                 ub_kind_set, torsion_kind_set, impr_kind_set, &
     777             :                                 opbend_kind_set, nelectron, nsgf, &
     778             :                                 molname_generated)
     779             : 
     780             :       TYPE(molecule_kind_type), INTENT(INOUT)            :: molecule_kind
     781             :       CHARACTER(LEN=*), INTENT(IN), OPTIONAL             :: name
     782             :       REAL(KIND=dp), OPTIONAL                            :: mass, charge
     783             :       INTEGER, INTENT(IN), OPTIONAL                      :: kind_number
     784             :       INTEGER, DIMENSION(:), OPTIONAL, POINTER           :: molecule_list
     785             :       TYPE(atom_type), DIMENSION(:), OPTIONAL, POINTER   :: atom_list
     786             :       INTEGER, INTENT(IN), OPTIONAL                      :: nbond
     787             :       TYPE(bond_type), DIMENSION(:), OPTIONAL, POINTER   :: bond_list
     788             :       INTEGER, INTENT(IN), OPTIONAL                      :: nbend
     789             :       TYPE(bend_type), DIMENSION(:), OPTIONAL, POINTER   :: bend_list
     790             :       INTEGER, INTENT(IN), OPTIONAL                      :: nub
     791             :       TYPE(ub_type), DIMENSION(:), OPTIONAL, POINTER     :: ub_list
     792             :       INTEGER, INTENT(IN), OPTIONAL                      :: nimpr
     793             :       TYPE(impr_type), DIMENSION(:), OPTIONAL, POINTER   :: impr_list
     794             :       INTEGER, INTENT(IN), OPTIONAL                      :: nopbend
     795             :       TYPE(opbend_type), DIMENSION(:), OPTIONAL, POINTER :: opbend_list
     796             :       INTEGER, INTENT(IN), OPTIONAL                      :: ntorsion
     797             :       TYPE(torsion_type), DIMENSION(:), OPTIONAL, &
     798             :          POINTER                                         :: torsion_list
     799             :       TYPE(fixd_constraint_type), DIMENSION(:), &
     800             :          OPTIONAL, POINTER                               :: fixd_list
     801             :       TYPE(colvar_counters), INTENT(IN), OPTIONAL        :: ncolv
     802             :       TYPE(colvar_constraint_type), DIMENSION(:), &
     803             :          OPTIONAL, POINTER                               :: colv_list
     804             :       INTEGER, INTENT(IN), OPTIONAL                      :: ng3x3
     805             :       TYPE(g3x3_constraint_type), DIMENSION(:), &
     806             :          OPTIONAL, POINTER                               :: g3x3_list
     807             :       INTEGER, INTENT(IN), OPTIONAL                      :: ng4x6, nfixd
     808             :       TYPE(g4x6_constraint_type), DIMENSION(:), &
     809             :          OPTIONAL, POINTER                               :: g4x6_list
     810             :       INTEGER, INTENT(IN), OPTIONAL                      :: nvsite
     811             :       TYPE(vsite_constraint_type), DIMENSION(:), &
     812             :          OPTIONAL, POINTER                               :: vsite_list
     813             :       INTEGER, INTENT(IN), OPTIONAL                      :: ng3x3_restraint, ng4x6_restraint, &
     814             :                                                             nfixd_restraint, nshell
     815             :       TYPE(shell_type), DIMENSION(:), OPTIONAL, POINTER  :: shell_list
     816             :       INTEGER, INTENT(IN), OPTIONAL                      :: nvsite_restraint
     817             :       TYPE(bond_kind_type), DIMENSION(:), OPTIONAL, &
     818             :          POINTER                                         :: bond_kind_set
     819             :       TYPE(bend_kind_type), DIMENSION(:), OPTIONAL, &
     820             :          POINTER                                         :: bend_kind_set
     821             :       TYPE(ub_kind_type), DIMENSION(:), OPTIONAL, &
     822             :          POINTER                                         :: ub_kind_set
     823             :       TYPE(torsion_kind_type), DIMENSION(:), OPTIONAL, &
     824             :          POINTER                                         :: torsion_kind_set
     825             :       TYPE(impr_kind_type), DIMENSION(:), OPTIONAL, &
     826             :          POINTER                                         :: impr_kind_set
     827             :       TYPE(opbend_kind_type), DIMENSION(:), OPTIONAL, &
     828             :          POINTER                                         :: opbend_kind_set
     829             :       INTEGER, INTENT(IN), OPTIONAL                      :: nelectron, nsgf
     830             :       LOGICAL, INTENT(IN), OPTIONAL                      :: molname_generated
     831             : 
     832             :       INTEGER                                            :: n
     833             : 
     834     2005801 :       IF (PRESENT(atom_list)) THEN
     835      264842 :          n = SIZE(atom_list)
     836      264842 :          molecule_kind%natom = n
     837      264842 :          molecule_kind%atom_list => atom_list
     838             :       END IF
     839     2005801 :       IF (PRESENT(molname_generated)) molecule_kind%molname_generated = molname_generated
     840     2005801 :       IF (PRESENT(name)) molecule_kind%name = name
     841     2005801 :       IF (PRESENT(mass)) molecule_kind%mass = mass
     842     2005801 :       IF (PRESENT(charge)) molecule_kind%charge = charge
     843     2005801 :       IF (PRESENT(kind_number)) molecule_kind%kind_number = kind_number
     844     2005801 :       IF (PRESENT(nbond)) molecule_kind%nbond = nbond
     845     2005801 :       IF (PRESENT(bond_list)) molecule_kind%bond_list => bond_list
     846     2005801 :       IF (PRESENT(nbend)) molecule_kind%nbend = nbend
     847     2005801 :       IF (PRESENT(nelectron)) molecule_kind%nelectron = nelectron
     848     2005801 :       IF (PRESENT(nsgf)) molecule_kind%nsgf = nsgf
     849     2005801 :       IF (PRESENT(bend_list)) molecule_kind%bend_list => bend_list
     850     2005801 :       IF (PRESENT(nub)) molecule_kind%nub = nub
     851     2005801 :       IF (PRESENT(ub_list)) molecule_kind%ub_list => ub_list
     852     2005801 :       IF (PRESENT(ntorsion)) molecule_kind%ntorsion = ntorsion
     853     2005801 :       IF (PRESENT(torsion_list)) molecule_kind%torsion_list => torsion_list
     854     2005801 :       IF (PRESENT(nimpr)) molecule_kind%nimpr = nimpr
     855     2005801 :       IF (PRESENT(impr_list)) molecule_kind%impr_list => impr_list
     856     2005801 :       IF (PRESENT(nopbend)) molecule_kind%nopbend = nopbend
     857     2005801 :       IF (PRESENT(opbend_list)) molecule_kind%opbend_list => opbend_list
     858     2005801 :       IF (PRESENT(ncolv)) molecule_kind%ncolv = ncolv
     859     2005801 :       IF (PRESENT(colv_list)) molecule_kind%colv_list => colv_list
     860     2005801 :       IF (PRESENT(ng3x3)) molecule_kind%ng3x3 = ng3x3
     861     2005801 :       IF (PRESENT(g3x3_list)) molecule_kind%g3x3_list => g3x3_list
     862     2005801 :       IF (PRESENT(ng4x6)) molecule_kind%ng4x6 = ng4x6
     863     2005801 :       IF (PRESENT(nvsite)) molecule_kind%nvsite = nvsite
     864     2005801 :       IF (PRESENT(nfixd)) molecule_kind%nfixd = nfixd
     865     2005801 :       IF (PRESENT(nfixd_restraint)) molecule_kind%nfixd_restraint = nfixd_restraint
     866     2005801 :       IF (PRESENT(ng3x3_restraint)) molecule_kind%ng3x3_restraint = ng3x3_restraint
     867     2005801 :       IF (PRESENT(ng4x6_restraint)) molecule_kind%ng4x6_restraint = ng4x6_restraint
     868     2005801 :       IF (PRESENT(nvsite_restraint)) molecule_kind%nvsite_restraint = nvsite_restraint
     869     2005801 :       IF (PRESENT(g4x6_list)) molecule_kind%g4x6_list => g4x6_list
     870     2005801 :       IF (PRESENT(vsite_list)) molecule_kind%vsite_list => vsite_list
     871     2005801 :       IF (PRESENT(fixd_list)) molecule_kind%fixd_list => fixd_list
     872     2005801 :       IF (PRESENT(bond_kind_set)) molecule_kind%bond_kind_set => bond_kind_set
     873     2005801 :       IF (PRESENT(bend_kind_set)) molecule_kind%bend_kind_set => bend_kind_set
     874     2005801 :       IF (PRESENT(ub_kind_set)) molecule_kind%ub_kind_set => ub_kind_set
     875     2005801 :       IF (PRESENT(torsion_kind_set)) molecule_kind%torsion_kind_set => torsion_kind_set
     876     2005801 :       IF (PRESENT(impr_kind_set)) molecule_kind%impr_kind_set => impr_kind_set
     877     2005801 :       IF (PRESENT(opbend_kind_set)) molecule_kind%opbend_kind_set => opbend_kind_set
     878     2005801 :       IF (PRESENT(nshell)) molecule_kind%nshell = nshell
     879     2005801 :       IF (PRESENT(shell_list)) molecule_kind%shell_list => shell_list
     880     2005801 :       IF (PRESENT(molecule_list)) THEN
     881      132421 :          n = SIZE(molecule_list)
     882      132421 :          molecule_kind%nmolecule = n
     883      132421 :          molecule_kind%molecule_list => molecule_list
     884             :       END IF
     885     2005801 :    END SUBROUTINE set_molecule_kind
     886             : 
     887             : ! **************************************************************************************************
     888             : !> \brief   Write a molecule kind data set to the output unit.
     889             : !> \param molecule_kind ...
     890             : !> \param output_unit ...
     891             : !> \date    24.09.2003
     892             : !> \author  MK
     893             : !> \version 1.0
     894             : ! **************************************************************************************************
     895        2228 :    SUBROUTINE write_molecule_kind(molecule_kind, output_unit)
     896             :       TYPE(molecule_kind_type), INTENT(IN)               :: molecule_kind
     897             :       INTEGER, INTENT(in)                                :: output_unit
     898             : 
     899             :       CHARACTER(LEN=default_string_length)               :: name
     900             :       INTEGER                                            :: iatom, imolecule, natom, nmolecule
     901             :       TYPE(atomic_kind_type), POINTER                    :: atomic_kind
     902             : 
     903        2228 :       IF (output_unit > 0) THEN
     904        2228 :          natom = SIZE(molecule_kind%atom_list)
     905        2228 :          nmolecule = SIZE(molecule_kind%molecule_list)
     906             : 
     907        2228 :          IF (natom == 1) THEN
     908         211 :             atomic_kind => molecule_kind%atom_list(1)%atomic_kind
     909         211 :             CALL get_atomic_kind(atomic_kind=atomic_kind, name=name)
     910             :             WRITE (UNIT=output_unit, FMT="(/,T2,I5,A,T36,A,A,T64,A)") &
     911         211 :                molecule_kind%kind_number, &
     912         211 :                ". Molecule kind: "//TRIM(molecule_kind%name), &
     913         422 :                "Atomic kind name:   ", TRIM(name)
     914             :             WRITE (UNIT=output_unit, FMT="(T9,A,L1,T55,A,T75,I6)") &
     915         211 :                "Automatic name: ", molecule_kind%molname_generated, &
     916         422 :                "Number of molecules:", nmolecule
     917             :          ELSE
     918             :             WRITE (UNIT=output_unit, FMT="(/,T2,I5,A,T50,A,T75,I6,/,T22,A)") &
     919        2017 :                molecule_kind%kind_number, &
     920        2017 :                ". Molecule kind: "//TRIM(molecule_kind%name), &
     921        2017 :                "Number of atoms:    ", natom, &
     922        4034 :                "Atom         Atomic kind name"
     923       16712 :             DO iatom = 1, natom
     924       14695 :                atomic_kind => molecule_kind%atom_list(iatom)%atomic_kind
     925       14695 :                CALL get_atomic_kind(atomic_kind=atomic_kind, name=name)
     926             :                WRITE (UNIT=output_unit, FMT="(T20,I6,(7X,A18))") &
     927       16712 :                   iatom, TRIM(name)
     928             :             END DO
     929             :             WRITE (UNIT=output_unit, FMT="(/,T9,A,L1)") &
     930        2017 :                "The name was automatically generated: ", &
     931        4034 :                molecule_kind%molname_generated
     932             :             WRITE (UNIT=output_unit, FMT="(T9,A,I6,/,T9,A,(T30,5I10))") &
     933        2017 :                "Number of molecules: ", nmolecule, "Molecule list:", &
     934       33824 :                (molecule_kind%molecule_list(imolecule), imolecule=1, nmolecule)
     935        2017 :             IF (molecule_kind%nbond > 0) &
     936             :                WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
     937        1741 :                "Number of bonds:       ", molecule_kind%nbond
     938        2017 :             IF (molecule_kind%nbend > 0) &
     939             :                WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
     940        1615 :                "Number of bends:       ", molecule_kind%nbend
     941        2017 :             IF (molecule_kind%nub > 0) &
     942             :                WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
     943         264 :                "Number of Urey-Bradley:", molecule_kind%nub
     944        2017 :             IF (molecule_kind%ntorsion > 0) &
     945             :                WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
     946        1122 :                "Number of torsions:    ", molecule_kind%ntorsion
     947        2017 :             IF (molecule_kind%nimpr > 0) &
     948             :                WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
     949         179 :                "Number of improper:    ", molecule_kind%nimpr
     950        2017 :             IF (molecule_kind%nopbend > 0) &
     951             :                WRITE (UNIT=output_unit, FMT="(1X,A30,I6)") &
     952           4 :                "Number of out opbends:    ", molecule_kind%nopbend
     953             :          END IF
     954             :       END IF
     955        2228 :    END SUBROUTINE write_molecule_kind
     956             : 
     957             : ! **************************************************************************************************
     958             : !> \brief   Write a moleculeatomic kind set data set to the output unit.
     959             : !> \param molecule_kind_set ...
     960             : !> \param subsys_section ...
     961             : !> \date    24.09.2003
     962             : !> \author  MK
     963             : !> \version 1.0
     964             : ! **************************************************************************************************
     965        9493 :    SUBROUTINE write_molecule_kind_set(molecule_kind_set, subsys_section)
     966             :       TYPE(molecule_kind_type), DIMENSION(:), INTENT(IN) :: molecule_kind_set
     967             :       TYPE(section_vals_type), INTENT(IN)                :: subsys_section
     968             : 
     969             :       CHARACTER(len=*), PARAMETER :: routineN = 'write_molecule_kind_set'
     970             : 
     971             :       INTEGER                                            :: handle, imolecule_kind, natom, nbend, &
     972             :                                                             nbond, nimpr, nmolecule, &
     973             :                                                             nmolecule_kind, nopbend, ntors, &
     974             :                                                             ntotal, nub, output_unit
     975             :       LOGICAL                                            :: all_single_atoms
     976             :       TYPE(cp_logger_type), POINTER                      :: logger
     977             : 
     978        9493 :       CALL timeset(routineN, handle)
     979             : 
     980        9493 :       NULLIFY (logger)
     981        9493 :       logger => cp_get_default_logger()
     982             :       output_unit = cp_print_key_unit_nr(logger, subsys_section, &
     983        9493 :                                          "PRINT%MOLECULES", extension=".Log")
     984        9493 :       IF (output_unit > 0) THEN
     985        2566 :          WRITE (UNIT=output_unit, FMT="(/,/,T2,A)") "MOLECULE KIND INFORMATION"
     986             : 
     987        2566 :          nmolecule_kind = SIZE(molecule_kind_set)
     988             : 
     989        2566 :          all_single_atoms = .TRUE.
     990       31515 :          DO imolecule_kind = 1, nmolecule_kind
     991       28949 :             natom = SIZE(molecule_kind_set(imolecule_kind)%atom_list)
     992       28949 :             nmolecule = SIZE(molecule_kind_set(imolecule_kind)%molecule_list)
     993       31515 :             IF (natom*nmolecule > 1) all_single_atoms = .FALSE.
     994             :          END DO
     995             : 
     996        2566 :          IF (all_single_atoms) THEN
     997             :             WRITE (UNIT=output_unit, FMT="(/,/,T2,A)") &
     998        1879 :                "All atoms are their own molecule, skipping detailed information"
     999             :          ELSE
    1000        2915 :             DO imolecule_kind = 1, nmolecule_kind
    1001        2915 :                CALL write_molecule_kind(molecule_kind_set(imolecule_kind), output_unit)
    1002             :             END DO
    1003             :          END IF
    1004             : 
    1005             :          CALL get_molecule_kind_set(molecule_kind_set=molecule_kind_set, &
    1006             :                                     nbond=nbond, &
    1007             :                                     nbend=nbend, &
    1008             :                                     nub=nub, &
    1009             :                                     ntorsion=ntors, &
    1010             :                                     nimpr=nimpr, &
    1011        2566 :                                     nopbend=nopbend)
    1012        2566 :          ntotal = nbond + nbend + nub + ntors + nimpr + nopbend
    1013        2566 :          IF (ntotal > 0) THEN
    1014             :             WRITE (UNIT=output_unit, FMT="(/,/,T2,A,T45,A30,I6)") &
    1015         585 :                "MOLECULE KIND SET INFORMATION", &
    1016        1170 :                "Total Number of bonds:       ", nbond
    1017             :             WRITE (UNIT=output_unit, FMT="(T45,A30,I6)") &
    1018         585 :                "Total Number of bends:       ", nbend
    1019             :             WRITE (UNIT=output_unit, FMT="(T45,A30,I6)") &
    1020         585 :                "Total Number of Urey-Bradley:", nub
    1021             :             WRITE (UNIT=output_unit, FMT="(T45,A30,I6)") &
    1022         585 :                "Total Number of torsions:    ", ntors
    1023             :             WRITE (UNIT=output_unit, FMT="(T45,A30,I6)") &
    1024         585 :                "Total Number of improper:    ", nimpr
    1025             :             WRITE (UNIT=output_unit, FMT="(T45,A30,I6)") &
    1026         585 :                "Total Number of opbends:    ", nopbend
    1027             :          END IF
    1028             :       END IF
    1029             :       CALL cp_print_key_finished_output(output_unit, logger, subsys_section, &
    1030        9493 :                                         "PRINT%MOLECULES")
    1031             : 
    1032        9493 :       CALL timestop(handle)
    1033             : 
    1034        9493 :    END SUBROUTINE write_molecule_kind_set
    1035             : 
    1036           0 : END MODULE molecule_kind_types

Generated by: LCOV version 1.15