1      subroutine hf3PEabc(Pabc, E, R0, IJK,
2     &    Nabc, La, Lb, Lc, Lg, Lg3, nint, lim_nint)
3c $Id$
4      implicit none
5#include "errquit.fh"
6c
7      integer La, Lb, Lc, Lg, Lg3, Nabc, nint, lim_nint
8      integer IJK(0:Lg,0:Lg,0:Lg)
9      double precision E(Nabc,3,0:Lg,0:La,0:Lb,0:Lc)
10      double precision R0(Nabc,Lg3)
11      double precision Pabc(*)
12      double precision pabc_int
13c
14      integer Nxyz(3)
15      integer La2, Lb2, Lc2, ica, icb, icc
16      integer Ia, Ja, Ka
17      integer Ib, Jb, Kb
18      integer Ic, Jc, Kc
19      integer Ip, Jp, Kp
20      integer np, mp
21c
22c Define the number of shell components on each center.
23
24      La2 = ((La+1)*(La+2))/2
25      Lb2 = ((Lb+1)*(Lb+2))/2
26      Lc2 = ((Lc+1)*(Lc+2))/2
27      nint = la2*lb2*lc2
28      if (nint.gt.lim_nint) then
29        write(6,*)' la = ',la
30        write(6,*)' lb = ',lb
31        write(6,*)' lc = ',lc
32        write(6,*)' nint     = ',nint
33        write(6,*)' lim_nint = ',lim_nint
34        call errquit('hf3PEabc: something hosed',911, INT_ERR)
35      endif
36
37c Loop over shell components.
38c
39      nint = 0
40c
41      do 00100 ica = 1,La2
42
43        call getNxyz(La,ica,Nxyz)
44        Ia = Nxyz(1)
45        Ja = Nxyz(2)
46        Ka = Nxyz(3)
47
48        do 00200 icb = 1,Lb2
49
50          call getNxyz(Lb,icb,Nxyz)
51          Ib = Nxyz(1)
52          Jb = Nxyz(2)
53          Kb = Nxyz(3)
54
55          do 00300 icc = 1, Lc2
56
57            call getNxyz(Lc,icc,Nxyz)
58            Ic = Nxyz(1)
59            Jc = Nxyz(2)
60            Kc = Nxyz(3)
61
62            if (nint.eq.lim_nint) then
63              write(6,*)' nint =   ',nint
64              write(6,*)' lim_nint ',lim_nint
65              call errquit
66     &            ('hf3PEabc: buffer too small ',911, MEM_ERR)
67            endif
68
69            nint = nint + 1
70
71
72            pabc_int = 0.0d00
73            do 00400 Ip = 0, (Ia + Ib + Ic)
74              do 00500 Jp = 0, (Ja + Jb +Jc)
75                do 00600 Kp = 0, (Ka + Kb +Kc)
76
77                  np = IJK(Ip, Jp, Kp)
78
79                  do 00700 mp = 1, Nabc
80                    pabc_int = pabc_int +
81     &                  E(mp,1,Ip,Ia,Ib,Ic) *
82     &                  E(mp,2,Jp,Ja,Jb,Jc) *
83     &                  E(mp,3,Kp,Ka,Kb,Kc) * R0(mp,np)
84
8500700             continue
8600600           continue
8700500         continue
8800400       continue
89            Pabc(nint) = pabc_int
9000300     continue
9100200   continue
9200100 continue
93c
94      end
95