LCOV - code coverage report
Current view: top level - src/subsys - molecule_types.F (source / functions) Hit Total Coverage
Test: CP2K Regtests (git:34ef472) Lines: 124 133 93.2 %
Date: 2024-04-26 08:30:29 Functions: 8 15 53.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 data structure for the molecule information.
      10             : !> \par History
      11             : !>      JGH (22.05.2004) add last_atom information
      12             : !>      Teodoro Laino [tlaino] 12.2008 - Preparing for VIRTUAL SITE constraints
      13             : !>                                       (patch by Marcel Baer)
      14             : !> \author MK (29.08.2003)
      15             : ! **************************************************************************************************
      16             : MODULE molecule_types
      17             : 
      18             :    USE colvar_types,                    ONLY: colvar_counters,&
      19             :                                               colvar_release,&
      20             :                                               colvar_type
      21             :    USE kinds,                           ONLY: dp
      22             :    USE molecule_kind_types,             ONLY: colvar_constraint_type,&
      23             :                                               fixd_constraint_type,&
      24             :                                               g3x3_constraint_type,&
      25             :                                               g4x6_constraint_type,&
      26             :                                               get_molecule_kind,&
      27             :                                               molecule_kind_type,&
      28             :                                               vsite_constraint_type
      29             : #include "../base/base_uses.f90"
      30             : 
      31             :    IMPLICIT NONE
      32             : 
      33             :    PRIVATE
      34             : 
      35             : ! *** Global parameters (in this module) ***
      36             : 
      37             :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'molecule_types'
      38             : 
      39             : ! *** Data types ***
      40             : ! **************************************************************************************************
      41             :    TYPE local_colvar_constraint_type
      42             :       LOGICAL                       :: init = .FALSE.
      43             :       TYPE(colvar_type), POINTER     :: colvar => NULL()
      44             :       TYPE(colvar_type), POINTER     :: colvar_old => NULL()
      45             :       REAL(KIND=dp)               :: lambda = 0.0_dp, sigma = 0.0_dp
      46             :    END TYPE local_colvar_constraint_type
      47             : 
      48             : ! **************************************************************************************************
      49             :    TYPE local_g3x3_constraint_type
      50             :       LOGICAL                       :: init = .FALSE.
      51             :       REAL(KIND=dp)               :: scale = 0.0_dp, scale_old = 0.0_dp, &
      52             :                                      imass1 = 0.0_dp, imass2 = 0.0_dp, imass3 = 0.0_dp
      53             :       REAL(KIND=dp), DIMENSION(3) :: fa = 0.0_dp, fb = 0.0_dp, fc = 0.0_dp, &
      54             :                                      f_roll1 = 0.0_dp, f_roll2 = 0.0_dp, f_roll3 = 0.0_dp, &
      55             :                                      ra_old = 0.0_dp, rb_old = 0.0_dp, rc_old = 0.0_dp, &
      56             :                                      va = 0.0_dp, vb = 0.0_dp, vc = 0.0_dp, &
      57             :                                      lambda = 0.0_dp, del_lambda = 0.0_dp, lambda_old = 0.0_dp, &
      58             :                                      r0_12 = 0.0_dp, r0_13 = 0.0_dp, r0_23 = 0.0_dp
      59             :       REAL(KIND=dp), DIMENSION(3, 3) :: amat = 0.0_dp
      60             :    END TYPE local_g3x3_constraint_type
      61             : 
      62             : ! **************************************************************************************************
      63             :    TYPE local_g4x6_constraint_type
      64             :       LOGICAL                       :: init = .FALSE.
      65             :       REAL(KIND=dp)               :: scale = 0.0_dp, scale_old = 0.0_dp, imass1 = 0.0_dp, &
      66             :                                      imass2 = 0.0_dp, imass3 = 0.0_dp, imass4 = 0.0_dp
      67             :       REAL(KIND=dp), DIMENSION(3) :: fa = 0.0_dp, fb = 0.0_dp, fc = 0.0_dp, fd = 0.0_dp, fe = 0.0_dp, ff = 0.0_dp, &
      68             :                                      f_roll1 = 0.0_dp, f_roll2 = 0.0_dp, f_roll3 = 0.0_dp, &
      69             :                                      f_roll4 = 0.0_dp, f_roll5 = 0.0_dp, f_roll6 = 0.0_dp, &
      70             :                                      ra_old = 0.0_dp, rb_old = 0.0_dp, rc_old = 0.0_dp, &
      71             :                                      rd_old = 0.0_dp, re_old = 0.0_dp, rf_old = 0.0_dp, &
      72             :                                      va = 0.0_dp, vb = 0.0_dp, vc = 0.0_dp, vd = 0.0_dp, ve = 0.0_dp, vf = 0.0_dp, &
      73             :                                      r0_12 = 0.0_dp, r0_13 = 0.0_dp, r0_14 = 0.0_dp, &
      74             :                                      r0_23 = 0.0_dp, r0_24 = 0.0_dp, r0_34 = 0.0_dp
      75             :       REAL(KIND=dp), DIMENSION(6)   :: lambda = 0.0_dp, del_lambda = 0.0_dp, lambda_old = 0.0_dp
      76             :       REAL(KIND=dp), DIMENSION(6, 6) :: amat = 0.0_dp
      77             :    END TYPE local_g4x6_constraint_type
      78             : 
      79             : ! **************************************************************************************************
      80             :    TYPE local_states_type
      81             :       INTEGER, DIMENSION(:), POINTER :: states => NULL() ! indices of Kohn-Sham states for molecule
      82             :       INTEGER                        :: nstates = 0 ! Kohn-Sham states for molecule
      83             :    END TYPE local_states_type
      84             : 
      85             : ! **************************************************************************************************
      86             :    TYPE local_constraint_type
      87             :       TYPE(local_colvar_constraint_type), DIMENSION(:), POINTER :: lcolv => NULL()
      88             :       TYPE(local_g3x3_constraint_type), DIMENSION(:), POINTER :: lg3x3 => NULL()
      89             :       TYPE(local_g4x6_constraint_type), DIMENSION(:), POINTER :: lg4x6 => NULL()
      90             :    END TYPE local_constraint_type
      91             : 
      92             : ! **************************************************************************************************
      93             :    TYPE global_constraint_type
      94             :       TYPE(colvar_counters)                    :: ncolv = colvar_counters()
      95             :       INTEGER                                  :: ntot = 0, nrestraint = 0
      96             :       INTEGER                                  :: ng3x3 = 0, ng3x3_restraint = 0
      97             :       INTEGER                                  :: ng4x6 = 0, ng4x6_restraint = 0
      98             :       INTEGER                                  :: nvsite = 0, nvsite_restraint = 0
      99             :       TYPE(fixd_constraint_type), DIMENSION(:), POINTER   :: fixd_list => NULL()
     100             :       TYPE(colvar_constraint_type), DIMENSION(:), POINTER :: colv_list => NULL()
     101             :       TYPE(g3x3_constraint_type), DIMENSION(:), POINTER   :: g3x3_list => NULL()
     102             :       TYPE(g4x6_constraint_type), DIMENSION(:), POINTER   :: g4x6_list => NULL()
     103             :       TYPE(vsite_constraint_type), DIMENSION(:), POINTER  :: vsite_list => NULL()
     104             :       TYPE(local_colvar_constraint_type), DIMENSION(:), POINTER :: lcolv => NULL()
     105             :       TYPE(local_g3x3_constraint_type), DIMENSION(:), POINTER :: lg3x3 => NULL()
     106             :       TYPE(local_g4x6_constraint_type), DIMENSION(:), POINTER :: lg4x6 => NULL()
     107             :    END TYPE global_constraint_type
     108             : 
     109             : ! **************************************************************************************************
     110             :    TYPE molecule_type
     111             :       TYPE(molecule_kind_type), POINTER    :: molecule_kind => NULL() ! pointer to molecule kind information
     112             :       TYPE(local_states_type), DIMENSION(:), POINTER   :: lmi => NULL() ! local (spin)-states information
     113             :       TYPE(local_constraint_type), POINTER :: lci => NULL() ! local molecule constraint info
     114             :       INTEGER                              :: first_atom = 0 ! global index of first atom in molecule
     115             :       INTEGER                              :: last_atom = 0 ! global index of last atom in molecule
     116             :       INTEGER                              :: first_shell = 0 ! global index of first shell atom in molecule
     117             :       INTEGER                              :: last_shell = 0 ! global index of last shell atom in molecule
     118             :    END TYPE molecule_type
     119             : 
     120             : ! *** Public data types ***
     121             : 
     122             :    PUBLIC :: local_colvar_constraint_type, &
     123             :              local_g3x3_constraint_type, &
     124             :              local_g4x6_constraint_type, &
     125             :              local_constraint_type, &
     126             :              local_states_type, &
     127             :              global_constraint_type, &
     128             :              molecule_type
     129             : 
     130             : ! *** Public subroutines ***
     131             : 
     132             :    PUBLIC :: deallocate_global_constraint, &
     133             :              allocate_molecule_set, &
     134             :              deallocate_molecule_set, &
     135             :              get_molecule, &
     136             :              set_molecule, &
     137             :              set_molecule_set, &
     138             :              molecule_of_atom, &
     139             :              get_molecule_set_info
     140             : 
     141             : CONTAINS
     142             : 
     143             : ! **************************************************************************************************
     144             : !> \brief   Deallocate a global constraint.
     145             : !> \param gci ...
     146             : !> \par History
     147             : !>      07.2003 created [fawzi]
     148             : !>      01.2014 moved from cp_subsys_release() into separate routine.
     149             : !> \author  Ole Schuett
     150             : ! **************************************************************************************************
     151        9374 :    SUBROUTINE deallocate_global_constraint(gci)
     152             :       TYPE(global_constraint_type), POINTER              :: gci
     153             : 
     154             :       INTEGER                                            :: i
     155             : 
     156        9374 :       IF (ASSOCIATED(gci)) THEN
     157             :          ! List of constraints
     158        8846 :          IF (ASSOCIATED(gci%colv_list)) THEN
     159         110 :             DO i = 1, SIZE(gci%colv_list)
     160         110 :                DEALLOCATE (gci%colv_list(i)%i_atoms)
     161             :             END DO
     162          44 :             DEALLOCATE (gci%colv_list)
     163             :          END IF
     164             : 
     165        8846 :          IF (ASSOCIATED(gci%g3x3_list)) &
     166           4 :             DEALLOCATE (gci%g3x3_list)
     167             : 
     168        8846 :          IF (ASSOCIATED(gci%g4x6_list)) &
     169           4 :             DEALLOCATE (gci%g4x6_list)
     170             : 
     171             :          ! Local information
     172        8846 :          IF (ASSOCIATED(gci%lcolv)) THEN
     173         110 :             DO i = 1, SIZE(gci%lcolv)
     174          66 :                CALL colvar_release(gci%lcolv(i)%colvar)
     175         110 :                CALL colvar_release(gci%lcolv(i)%colvar_old)
     176             :             END DO
     177          44 :             DEALLOCATE (gci%lcolv)
     178             :          END IF
     179             : 
     180        8846 :          IF (ASSOCIATED(gci%lg3x3)) &
     181           4 :             DEALLOCATE (gci%lg3x3)
     182             : 
     183        8846 :          IF (ASSOCIATED(gci%lg4x6)) &
     184           4 :             DEALLOCATE (gci%lg4x6)
     185             : 
     186        8846 :          IF (ASSOCIATED(gci%fixd_list)) &
     187           2 :             DEALLOCATE (gci%fixd_list)
     188             : 
     189        8846 :          DEALLOCATE (gci)
     190             :       END IF
     191        9374 :    END SUBROUTINE deallocate_global_constraint
     192             : 
     193             : ! **************************************************************************************************
     194             : !> \brief   Allocate a molecule set.
     195             : !> \param molecule_set ...
     196             : !> \param nmolecule ...
     197             : !> \date    29.08.2003
     198             : !> \author  MK
     199             : !> \version 1.0
     200             : ! **************************************************************************************************
     201        9374 :    SUBROUTINE allocate_molecule_set(molecule_set, nmolecule)
     202             :       TYPE(molecule_type), DIMENSION(:), POINTER         :: molecule_set
     203             :       INTEGER, INTENT(IN)                                :: nmolecule
     204             : 
     205        9374 :       IF (ASSOCIATED(molecule_set)) CALL deallocate_molecule_set(molecule_set)
     206             : 
     207      328454 :       ALLOCATE (molecule_set(nmolecule))
     208             : 
     209        9374 :    END SUBROUTINE allocate_molecule_set
     210             : 
     211             : ! **************************************************************************************************
     212             : !> \brief   Deallocate a molecule set.
     213             : !> \param molecule_set ...
     214             : !> \date    29.08.2003
     215             : !> \author  MK
     216             : !> \version 1.0
     217             : ! **************************************************************************************************
     218        9374 :    SUBROUTINE deallocate_molecule_set(molecule_set)
     219             :       TYPE(molecule_type), DIMENSION(:), POINTER         :: molecule_set
     220             : 
     221             :       INTEGER                                            :: imolecule, j
     222             : 
     223        9374 :       IF (ASSOCIATED(molecule_set)) THEN
     224             : 
     225      309706 :          DO imolecule = 1, SIZE(molecule_set)
     226      300332 :             IF (ASSOCIATED(molecule_set(imolecule)%lmi)) THEN
     227          49 :                DO j = 1, SIZE(molecule_set(imolecule)%lmi)
     228          49 :                   IF (ASSOCIATED(molecule_set(imolecule)%lmi(j)%states)) THEN
     229          28 :                      DEALLOCATE (molecule_set(imolecule)%lmi(j)%states)
     230             :                   END IF
     231             :                END DO
     232          21 :                DEALLOCATE (molecule_set(imolecule)%lmi)
     233             :             END IF
     234      309706 :             IF (ASSOCIATED(molecule_set(imolecule)%lci)) THEN
     235       43692 :                IF (ASSOCIATED(molecule_set(imolecule)%lci%lcolv)) THEN
     236        4336 :                   DO j = 1, SIZE(molecule_set(imolecule)%lci%lcolv)
     237        2228 :                      CALL colvar_release(molecule_set(imolecule)%lci%lcolv(j)%colvar)
     238        4336 :                      CALL colvar_release(molecule_set(imolecule)%lci%lcolv(j)%colvar_old)
     239             :                   END DO
     240        2108 :                   DEALLOCATE (molecule_set(imolecule)%lci%lcolv)
     241             :                END IF
     242       43692 :                IF (ASSOCIATED(molecule_set(imolecule)%lci%lg3x3)) THEN
     243       36354 :                   DEALLOCATE (molecule_set(imolecule)%lci%lg3x3)
     244             :                END IF
     245       43692 :                IF (ASSOCIATED(molecule_set(imolecule)%lci%lg4x6)) THEN
     246         650 :                   DEALLOCATE (molecule_set(imolecule)%lci%lg4x6)
     247             :                END IF
     248       43692 :                DEALLOCATE (molecule_set(imolecule)%lci)
     249             :             END IF
     250             :          END DO
     251        9374 :          DEALLOCATE (molecule_set)
     252             : 
     253             :       END IF
     254        9374 :       NULLIFY (molecule_set)
     255             : 
     256        9374 :    END SUBROUTINE deallocate_molecule_set
     257             : 
     258             : ! **************************************************************************************************
     259             : !> \brief   Get components from a molecule data set.
     260             : !> \param molecule ...
     261             : !> \param molecule_kind ...
     262             : !> \param lmi ...
     263             : !> \param lci ...
     264             : !> \param lg3x3 ...
     265             : !> \param lg4x6 ...
     266             : !> \param lcolv ...
     267             : !> \param first_atom ...
     268             : !> \param last_atom ...
     269             : !> \param first_shell ...
     270             : !> \param last_shell ...
     271             : !> \date    29.08.2003
     272             : !> \author  MK
     273             : !> \version 1.0
     274             : ! **************************************************************************************************
     275     8747513 :    SUBROUTINE get_molecule(molecule, molecule_kind, lmi, lci, lg3x3, lg4x6, lcolv, &
     276             :                            first_atom, last_atom, first_shell, last_shell)
     277             : 
     278             :       TYPE(molecule_type), INTENT(IN)                    :: molecule
     279             :       TYPE(molecule_kind_type), OPTIONAL, POINTER        :: molecule_kind
     280             :       TYPE(local_states_type), DIMENSION(:), OPTIONAL, &
     281             :          POINTER                                         :: lmi
     282             :       TYPE(local_constraint_type), OPTIONAL, POINTER     :: lci
     283             :       TYPE(local_g3x3_constraint_type), OPTIONAL, &
     284             :          POINTER                                         :: lg3x3(:)
     285             :       TYPE(local_g4x6_constraint_type), OPTIONAL, &
     286             :          POINTER                                         :: lg4x6(:)
     287             :       TYPE(local_colvar_constraint_type), DIMENSION(:), &
     288             :          OPTIONAL, POINTER                               :: lcolv
     289             :       INTEGER, OPTIONAL                                  :: first_atom, last_atom, first_shell, &
     290             :                                                             last_shell
     291             : 
     292     8747513 :       IF (PRESENT(first_atom)) first_atom = molecule%first_atom
     293     8747513 :       IF (PRESENT(last_atom)) last_atom = molecule%last_atom
     294     8747513 :       IF (PRESENT(first_shell)) first_shell = molecule%first_shell
     295     8747513 :       IF (PRESENT(last_shell)) last_shell = molecule%last_shell
     296     8747513 :       IF (PRESENT(molecule_kind)) molecule_kind => molecule%molecule_kind
     297     8747513 :       IF (PRESENT(lmi)) lmi => molecule%lmi
     298     8747513 :       IF (PRESENT(lci)) lci => molecule%lci
     299     8747513 :       IF (PRESENT(lcolv)) THEN
     300      928471 :          IF (ASSOCIATED(molecule%lci)) THEN
     301      928471 :             lcolv => molecule%lci%lcolv
     302             :          ELSE
     303           0 :             CPABORT("The pointer lci is not associated")
     304             :          END IF
     305             :       END IF
     306     8747513 :       IF (PRESENT(lg3x3)) THEN
     307     1530904 :          IF (ASSOCIATED(molecule%lci)) THEN
     308     1530904 :             lg3x3 => molecule%lci%lg3x3
     309             :          ELSE
     310           0 :             CPABORT("The pointer lci is not associated")
     311             :          END IF
     312             :       END IF
     313     8747513 :       IF (PRESENT(lg4x6)) THEN
     314      885788 :          IF (ASSOCIATED(molecule%lci)) THEN
     315      885788 :             lg4x6 => molecule%lci%lg4x6
     316             :          ELSE
     317           0 :             CPABORT("The pointer lci is not associated")
     318             :          END IF
     319             :       END IF
     320             : 
     321     8747513 :    END SUBROUTINE get_molecule
     322             : 
     323             : ! **************************************************************************************************
     324             : !> \brief   Set a molecule data set.
     325             : !> \param molecule ...
     326             : !> \param molecule_kind ...
     327             : !> \param lmi ...
     328             : !> \param lci ...
     329             : !> \param lcolv ...
     330             : !> \param lg3x3 ...
     331             : !> \param lg4x6 ...
     332             : !> \date    29.08.2003
     333             : !> \author  MK
     334             : !> \version 1.0
     335             : ! **************************************************************************************************
     336      936444 :    SUBROUTINE set_molecule(molecule, molecule_kind, lmi, lci, lcolv, lg3x3, lg4x6)
     337             :       TYPE(molecule_type), INTENT(INOUT)                 :: molecule
     338             :       TYPE(molecule_kind_type), OPTIONAL, POINTER        :: molecule_kind
     339             :       TYPE(local_states_type), DIMENSION(:), OPTIONAL, &
     340             :          POINTER                                         :: lmi
     341             :       TYPE(local_constraint_type), OPTIONAL, POINTER     :: lci
     342             :       TYPE(local_colvar_constraint_type), DIMENSION(:), &
     343             :          OPTIONAL, POINTER                               :: lcolv
     344             :       TYPE(local_g3x3_constraint_type), OPTIONAL, &
     345             :          POINTER                                         :: lg3x3(:)
     346             :       TYPE(local_g4x6_constraint_type), OPTIONAL, &
     347             :          POINTER                                         :: lg4x6(:)
     348             : 
     349      936444 :       IF (PRESENT(molecule_kind)) molecule%molecule_kind => molecule_kind
     350      936444 :       IF (PRESENT(lmi)) molecule%lmi => lmi
     351      936444 :       IF (PRESENT(lci)) molecule%lci => lci
     352      936444 :       IF (PRESENT(lcolv)) THEN
     353        2108 :          IF (ASSOCIATED(molecule%lci)) THEN
     354        2108 :             molecule%lci%lcolv => lcolv
     355             :          ELSE
     356           0 :             CPABORT("The pointer lci is not associated")
     357             :          END IF
     358             :       END IF
     359      936444 :       IF (PRESENT(lg3x3)) THEN
     360       36354 :          IF (ASSOCIATED(molecule%lci)) THEN
     361       36354 :             molecule%lci%lg3x3 => lg3x3
     362             :          ELSE
     363           0 :             CPABORT("The pointer lci is not associated")
     364             :          END IF
     365             :       END IF
     366      936444 :       IF (PRESENT(lg4x6)) THEN
     367         650 :          IF (ASSOCIATED(molecule%lci)) THEN
     368         650 :             molecule%lci%lg4x6 => lg4x6
     369             :          ELSE
     370           0 :             CPABORT("The pointer lci is not associated")
     371             :          END IF
     372             :       END IF
     373             : 
     374      936444 :    END SUBROUTINE set_molecule
     375             : 
     376             : ! **************************************************************************************************
     377             : !> \brief   Set a molecule data set.
     378             : !> \param molecule_set ...
     379             : !> \param first_atom ...
     380             : !> \param last_atom ...
     381             : !> \date    29.08.2003
     382             : !> \author  MK
     383             : !> \version 1.0
     384             : ! **************************************************************************************************
     385        9374 :    SUBROUTINE set_molecule_set(molecule_set, first_atom, last_atom)
     386             :       TYPE(molecule_type), DIMENSION(:), INTENT(INOUT)   :: molecule_set
     387             :       INTEGER, DIMENSION(:), INTENT(IN), OPTIONAL        :: first_atom, last_atom
     388             : 
     389             :       INTEGER                                            :: imolecule
     390             : 
     391        9374 :       IF (PRESENT(first_atom)) THEN
     392        9374 :          IF (SIZE(first_atom) /= SIZE(molecule_set)) THEN
     393             :             CALL cp_abort(__LOCATION__, &
     394             :                           "The sizes of first_atom and molecule_set "// &
     395           0 :                           "are different")
     396             :          END IF
     397             : 
     398      309706 :          DO imolecule = 1, SIZE(molecule_set)
     399      309706 :             molecule_set(imolecule)%first_atom = first_atom(imolecule)
     400             :          END DO
     401             :       END IF
     402             : 
     403        9374 :       IF (PRESENT(last_atom)) THEN
     404        9374 :          IF (SIZE(last_atom) /= SIZE(molecule_set)) THEN
     405             :             CALL cp_abort(__LOCATION__, &
     406             :                           "The sizes of last_atom and molecule_set "// &
     407           0 :                           "are different")
     408             :          END IF
     409             : 
     410      309706 :          DO imolecule = 1, SIZE(molecule_set)
     411      309706 :             molecule_set(imolecule)%last_atom = last_atom(imolecule)
     412             :          END DO
     413             :       END IF
     414             : 
     415        9374 :    END SUBROUTINE set_molecule_set
     416             : 
     417             : ! **************************************************************************************************
     418             : !> \brief   finds for each atom the molecule it belongs to
     419             : !> \param molecule_set ...
     420             : !> \param atom_to_mol ...
     421             : ! **************************************************************************************************
     422         458 :    SUBROUTINE molecule_of_atom(molecule_set, atom_to_mol)
     423             :       TYPE(molecule_type), DIMENSION(:), INTENT(IN)      :: molecule_set
     424             :       INTEGER, DIMENSION(:), INTENT(OUT)                 :: atom_to_mol
     425             : 
     426             :       INTEGER                                            :: first_atom, iatom, imol, last_atom
     427             : 
     428        3092 :       DO imol = 1, SIZE(molecule_set)
     429        2634 :          CALL get_molecule(molecule=molecule_set(imol), first_atom=first_atom, last_atom=last_atom)
     430        9074 :          DO iatom = first_atom, last_atom
     431        8616 :             atom_to_mol(iatom) = imol
     432             :          END DO ! iatom
     433             :       END DO ! imol
     434             : 
     435         458 :    END SUBROUTINE molecule_of_atom
     436             : 
     437             : ! **************************************************************************************************
     438             : !> \brief returns information about molecules in the set.
     439             : !> \param molecule_set ...
     440             : !> \param atom_to_mol ...
     441             : !> \param mol_to_first_atom ...
     442             : !> \param mol_to_last_atom ...
     443             : !> \param mol_to_nelectrons ...
     444             : !> \param mol_to_nbasis ...
     445             : !> \param mol_to_charge ...
     446             : !> \param mol_to_multiplicity ...
     447             : !> \par History
     448             : !>       2011.06 created [Rustam Z Khaliullin]
     449             : !> \author Rustam Z Khaliullin
     450             : ! **************************************************************************************************
     451         726 :    SUBROUTINE get_molecule_set_info(molecule_set, atom_to_mol, mol_to_first_atom, &
     452         726 :                                     mol_to_last_atom, mol_to_nelectrons, mol_to_nbasis, mol_to_charge, &
     453         242 :                                     mol_to_multiplicity)
     454             : 
     455             :       TYPE(molecule_type), DIMENSION(:), INTENT(IN)      :: molecule_set
     456             :       INTEGER, DIMENSION(:), INTENT(OUT), OPTIONAL :: atom_to_mol, mol_to_first_atom, &
     457             :          mol_to_last_atom, mol_to_nelectrons, mol_to_nbasis, mol_to_charge, mol_to_multiplicity
     458             : 
     459             :       INTEGER                                            :: first_atom, iatom, imol, last_atom, &
     460             :                                                             nbasis, nelec
     461             :       REAL(KIND=dp)                                      :: charge
     462             :       TYPE(molecule_kind_type), POINTER                  :: imol_kind
     463             : 
     464        1894 :       DO imol = 1, SIZE(molecule_set)
     465             : 
     466             :          CALL get_molecule(molecule=molecule_set(imol), molecule_kind=imol_kind, &
     467        1652 :                            first_atom=first_atom, last_atom=last_atom)
     468             : 
     469        1652 :          IF (PRESENT(mol_to_nelectrons)) THEN
     470         810 :             CALL get_molecule_kind(imol_kind, nelectron=nelec)
     471         810 :             mol_to_nelectrons(imol) = nelec
     472             :          END IF
     473             : 
     474        1652 :          IF (PRESENT(mol_to_multiplicity)) THEN
     475             :             ! RZK-warning: At the moment we can only get the total number
     476             :             !  of electrons (alpha+beta) and we do not have a way to get the multiplicity of mols.
     477             :             !  Therefore, the best we can do is to assume the singlet state for even number of electrons
     478             :             !  and doublet state for odd number of electrons (assume ne_alpha > ne_beta).
     479             :             !  The best way to implement a correct multiplicity subroutine in the future is to get
     480             :             !  the number of alpha and beta e- for each atom from init_atom_electronic_state. This way (as opposed to
     481             :             !  reading the multiplicities from file) the number of occupied and virtual orbitals
     482             :             !  will be consistent with atomic guess. A guess with broken symmetry will be easy to
     483             :             !  implement as well.
     484         842 :             CALL get_molecule_kind(imol_kind, nelectron=nelec)
     485         842 :             IF (MOD(nelec, 2) .EQ. 0) THEN
     486         838 :                mol_to_multiplicity(imol) = 1
     487             :             ELSE
     488           4 :                mol_to_multiplicity(imol) = 2
     489             :             END IF
     490             :          END IF
     491             : 
     492        1652 :          IF (PRESENT(mol_to_charge)) THEN
     493         842 :             CALL get_molecule_kind(imol_kind, charge=charge)
     494         842 :             mol_to_charge(imol) = NINT(charge)
     495             :          END IF
     496             : 
     497        1652 :          IF (PRESENT(mol_to_nbasis)) THEN
     498         810 :             CALL get_molecule_kind(imol_kind, nsgf=nbasis)
     499         810 :             mol_to_nbasis(imol) = nbasis
     500             :          END IF
     501             : 
     502        1652 :          IF (PRESENT(mol_to_first_atom)) THEN
     503        1652 :             mol_to_first_atom(imol) = first_atom
     504             :          END IF
     505             : 
     506        1652 :          IF (PRESENT(mol_to_last_atom)) THEN
     507        1652 :             mol_to_last_atom(imol) = last_atom
     508             :          END IF
     509             : 
     510        1894 :          IF (PRESENT(atom_to_mol)) THEN
     511        2766 :             DO iatom = first_atom, last_atom
     512        2766 :                atom_to_mol(iatom) = imol
     513             :             END DO ! iatom
     514             :          END IF
     515             : 
     516             :       END DO ! imol
     517             : 
     518         242 :    END SUBROUTINE get_molecule_set_info
     519             : 
     520           0 : END MODULE molecule_types

Generated by: LCOV version 1.15