LCOV - code coverage report
Current view: top level - src/subsys - molecule_kind_types.F (source / functions) Hit Total Coverage
Test: CP2K Regtests (git:20fe009) Lines: 352 364 96.7 %
Date: 2022-07-05 19:56:53 Functions: 8 24 33.3 %

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

Generated by: LCOV version 1.15