Statistics
| Revision:

root / src / CalcCnct.f90 @ 10

History | View | Annotate | Download (3.9 kB)

1
!================================================================
2
! Calculation of the connectivity.
3
!
4
! Input:
5
! Na: (INTEGER) Number of atoms in the system
6
! Atome: (INTEGER) Number of mass of the atoms.
7
!                  Needed to find their covalent raii.
8
!  x(Nat), y(nat), z(nat) : (REAL) Cartesian coordinates of the system.
9
!  LIAISONS: (INTEGER Na,0:NMaxL) Array containing the connectivity.
10
!            That is Liaisons(i,:) contains the atoms linked to i.
11
!            Liaisons(i,0) contains the number of atoms linked to i.
12
! r_cov(0:Max_Z) (REAL) Covalent radii
13
! Fact (REAL) Factor used to determined if two atoms are linked :
14
!            i and j are linked if (dist(i,j)*fact <= rcov(i)+rcov(j)
15
!
16
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
17
! Taken from Path_Module
18
! NMaxL (INTEGER): Maximum number of bonds for one atome
19
! Max_Z (INTEGER): Last atom for which I have the atomic of mass and covalent radii.
20
! Nom(Nat) (STRING): Name of the atoms.
21
! Prog (SCHARS): Name of the program used to calculate energies and gradients
22
!
23
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
24
! v1.0 (c) PFL
25
! v1.1 (c) PFL 11.2007
26
! A test has been added for periodic system, ie the linked test
27
! is done between i in the central cell, and j in all surrounding cells.
28
!
29
!================================================================
30

    
31
SUBROUTINE CalcCnct(na,atome,x,y,z,LIAISONS,r_cov,fact)
32

    
33
  use Path_module, only : NMaxL, max_Z,Nom,Prog, KINT, KREAL, FPBC, &
34
       kaBeg,kaEnd,kbBeg,kbEnd,kcBeg,kcEnd
35
  
36
  IMPLICIT NONE
37

    
38
  integer(KINT) :: na, atome(Na)
39
  real(KREAL) ::  x(Na),y(Na),z(Na),fact
40
  real(KREAL) :: vx,vy,vz,dist
41
  REAL(KREAL) :: r_cov(0:Max_Z)
42
  INTEGER(KINT) :: LIAISONS(Na,0:NMaxL),Nbli,Nblj
43

    
44
  ! For periodic systems
45
  INTEGER(KINT) :: Ka,Kb,Kc
46
  LOGICAL   :: Bound
47

    
48
  ! Internals
49
  logical(KREAL) ::debug
50
  REAL(KREAL) :: DistTh
51
  INTEGER(KINT) :: I,j,iat
52

    
53
  INTERFACE
54
     function valid(string) result (isValid)
55
       CHARACTER(*), intent(in) :: string
56
       logical                  :: isValid
57
     END function VALID
58
  END INTERFACE
59

    
60
  debug=valid('calccnct')
61
  if (debug)  Call Header (" Entering CalcCnct ")
62

    
63
  DO i=1,na
64
     Liaisons(i,0)=0
65
  END DO
66

    
67
  if (debug) THEN
68
     WRITE(*,*) 'CalcCnct : covalent radii used'
69
     DO iat=1,na
70
        i=atome(iat)
71
        WRITE(*,*) Nom(I),I,r_cov(i)*fact
72
     END DO
73
     WRITE(*,*) 'Coordinates'
74
     DO iat=1,na
75
        i=atome(iat)
76
        WRITE(*,*) Nom(I),x(iat),y(iat),z(iat)
77
     END DO
78
  END IF
79

    
80
  IF (FPBC) THEN
81
     DO i=1,na
82
        NbLi=LIAISONS(i,0)
83
        DO j=i+1,na
84
           Bound=.FALSE.
85
           DO ka=kaBeg,kaEnd
86
              DO Kb=kbBeg,kbEnd
87
                 DO Kc=kcBeg,kcEnd
88
                    CALL VectorPer(j,i,ka,kb,kc,x,y,z,vx,vy,vz,dist)
89
                    dist=dist/fact
90
                    distth=(r_cov(atome(i))+r_cov(atome(j)))/100.
91
!                    if (debug) WRITE(*,*) atome(i),atome(j),dist,distth
92
                    if (dist.le.distth) Bound=.TRUE.
93
                 END DO
94
              END DO
95
           END DO
96
           IF (Bound) THEN
97
              if (debug) WRITE(*,*) "Adding a bond between:",i,j
98
              NbLi=NbLi+1
99
              LIAISONS(i,NbLi)=j;
100
              NBlj=LIAISONS(j,0)+1
101
              LIAISONS(j,Nblj)=i;
102
              LIAISONS(j,0)=Nblj
103
           END IF
104
        END DO
105
        LIAISONS(i,0)=Nbli
106
     END DO
107
  ELSE
108
     DO i=1,na
109
        NbLi=LIAISONS(i,0)
110
        DO j=i+1,na
111
           CALL vecteur(j,i,x,y,z,vx,vy,vz,dist)
112
           dist=dist/fact
113
           distth=(r_cov(atome(i))+r_cov(atome(j)))/100.
114
!           if (debug) WRITE(*,*) atome(i),atome(j),dist,distth
115
           if (dist.le.distth) THEN
116
              NbLi=NbLi+1
117
              LIAISONS(i,NbLi)=j;
118
              NBlj=LIAISONS(j,0)+1
119
              LIAISONS(j,Nblj)=i;
120
              LIAISONS(j,0)=Nblj
121
           END IF
122
        END DO
123
        LIAISONS(i,0)=Nbli
124
     END DO
125
  END IF
126
  if (debug) Call Header (" CalcCnct over ")
127
END SUBROUTINE CalcCnct