1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
|
SUBROUTINE DIPIND (DIPVEC)
C...............................................................
C MODIFICATION OF DIPOLE SUBROUTINE FOR USE IN THE CALCULATION
C OF THE INDUCED DIPOLES FOR POLARIZABILITIES.
C...............................................................
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
INCLUDE 'SIZES'
COMMON /CORE / CORE(107)
COMMON /DENSTY/ P(MPACK),PA(MPACK),PB(MPACK)
COMMON /MOLMEC/ HTYPE(4),NHCO(4,20),NNHCO,ITYPE
COMMON /GEOM / GEO(3,NUMATM), XCOORD(3,NUMATM)
COMMON /MOLKST/ NUMAT,NAT(NUMATM),NFIRST(NUMATM),NMIDLE(NUMATM),
1 NLAST(NUMATM),NORBS,NELECS,NALPHA,NBETA,
2 NCLOSE,NOPEN,NDUMY,FRACT
COMMON /NUMCAL/ NUMCAL
COMMON /KEYWRD/ KEYWRD
COMMON /ISTOPE/ AMS(107)
COMMON /MULTIP/ DD(107), QQ(107), AM(107), AD(107), AQ(107)
DIMENSION Q(MAXORB),Q2(MAXORB),DIPVEC(3),CENTER(3),
1 COORD(3,NUMATM)
CHARACTER*241 KEYWRD
C
C***********************************************************************
C DIPOLE CALCULATES DIPOLE MOMENTS
C
C ON INPUT P = DENSITY MATRIX
C Q = TOTAL ATOMIC CHARGES, (NUCLEAR + ELECTRONIC)
C NUMAT = NUMBER OF ATOMS IN MOLECULE
C NAT = ATOMIC NUMBERS OF ATOMS
C NFIRST= START OF ATOM ORBITAL COUNTERS
C COORD = COORDINATES OF ATOMS
C
C OUTPUT DIPOLE = DIPOLE MOMENT
C***********************************************************************
C
C IN THE ZDO APPROXIMATION, ONLY TWO TERMS ARE RETAINED IN THE
C CALCULATION OF DIPOLE MOMENTS.
C 1. THE POINT CHARGE TERM (INDEPENDENT OF PARAMETERIZATION).
C 2. THE ONE-CENTER HYBRIDIZATION TERM, WHICH ARISES FROM MATRIX
C ELEMENTS OF THE FORM <NS/R/NP>. THIS TERM IS A FUNCTION OF
C THE SLATER EXPONENTS (ZS,ZP) AND IS THUS DEPENDENT ON PARAMETER-
C IZATION. THE HYBRIDIZATION FACTORS (HYF(I)) USED IN THIS SUB-
C ROUTINE ARE CALCULATED FROM THE FOLLOWING FORMULAE.
C FOR SECOND ROW ELEMENTS <2S/R/2P>
C HYF(I)= 469.56193322*(SQRT(((ZS(I)**5)*(ZP(I)**5)))/
C ((ZS(I) + ZP(I))**6))
C FOR THIRD ROW ELEMENTS <3S/R/3P>
C HYF(I)=2629.107682607*(SQRT(((ZS(I)**7)*(ZP(I)**7)))/
C ((ZS(I) + ZP(I))**8))
C FOR FOURTH ROW ELEMENTS AND UP :
C HYF(I)=2*(2.10716)*DD(I)
C WHERE DD(I) IS THE CHARGE SEPARATION IN ATOMIC UNITS
C
C
C REFERENCES:
C J.A.POPLE & D.L.BEVERIDGE: APPROXIMATE M.O. THEORY
C S.P.MCGLYNN, ET AL: APPLIED QUANTUM CHEMISTRY
C
DIMENSION DIP(4,3)
DIMENSION HYF(107,2)
SAVE ICALCN, HYF, WTMOL, CHARGD
LOGICAL CHARGD
DATA HYF(1,1) / 0.0D00 /
DATA HYF(1,2) /0.0D0 /
DATA HYF(5,2) /6.520587D0/
DATA HYF(6,2) /4.253676D0/
DATA HYF(7,2) /2.947501D0/
DATA HYF(8,2) /2.139793D0/
DATA HYF(9,2) /2.2210719D0/
DATA HYF(14,2)/6.663059D0/
DATA HYF(15,2)/5.657623D0/
DATA HYF(16,2)/6.345552D0/
DATA HYF(17,2)/2.522964D0/
DATA ICALCN/0/
C
C SETUP FOR DIPOLE CALCULATION
C
CALL CHRGE (P,Q2)
DO 10 I = 1,NUMAT
Q(I) = CORE(NAT(I)) - Q2(I)
10 CONTINUE
CALL GMETRY (GEO,COORD)
C
IF (ICALCN.NE.NUMCAL) THEN
DO 20 I=2,107
20 HYF(I,1)= 5.0832D0*DD(I)
WTMOL=0.D0
SUM=0.D0
DO 30 I=1,NUMAT
WTMOL=WTMOL+AMS(NAT(I))
30 SUM=SUM+Q(I)
CHARGD=(ABS(SUM).GT.0.5D0)
ICALCN=NUMCAL
KTYPE=1
IF(ITYPE.EQ.4)KTYPE=2
ENDIF
IF(CHARGD)THEN
C
C NEED TO RESET ION'S POSITION SO THAT THE CENTER OF MASS IS AT THE
C ORIGIN.
C
C$DOIT ASIS
DO 40 I=1,3
40 CENTER(I)=0.D0
DO 50 I=1,3
C$DOIT VBEST
DO 50 J=1,NUMAT
50 CENTER(I)=CENTER(I)+AMS(NAT(J))*COORD(I,J)
C$DOIT ASIS
DO 60 I=1,3
60 CENTER(I)=CENTER(I)/WTMOL
DO 70 I=1,3
C$DOIT VBEST
DO 70 J=1,NUMAT
70 COORD(I,J)=COORD(I,J)-CENTER(I)
ENDIF
C$DOIT ASIS
DO 80 I=1,4
C$DOIT ASIS
DO 80 J=1,3
80 DIP(I,J)=0.0D00
C$DOIT ASIS
DO 100 I=1,NUMAT
NI=NAT(I)
IA=NFIRST(I)
L=NLAST(I)-IA
C$DOIT ASIS
DO 90 J=1,L
K=((IA+J)*(IA+J-1))/2+IA
90 DIP(J,2)=DIP(J,2)-HYF(NI,KTYPE)*P(K)
C$DOIT ASIS
DO 100 J=1,3
100 DIP(J,1)=DIP(J,1)+4.803D00*Q(I)*COORD(J,I)
C$DOIT ASIS
DO 110 J=1,3
110 DIP(J,3)=DIP(J,2)+DIP(J,1)
C$DOIT ASIS
DO 120 J=1,3
120 DIP(4,J)=SQRT(DIP(1,J)**2+DIP(2,J)**2+DIP(3,J)**2)
DIPVEC(1)= -DIP(1,3)
DIPVEC(2)= -DIP(2,3)
DIPVEC(3)= -DIP(3,3)
C WRITE (6,60) ((DIP(I,J),I=1,4),J=1,3)
C 60 FORMAT (3(4F10.3))
RETURN
C
END
|
