1C> \ingroup selci
2C> @{
3C>
4      subroutine selci_makhdb(hij,ind,icase,ns,indbar,iocc,
5     $     w1,w2,work1,work2,g,int12,int34,numf,numf2)
6*
7* $Id$
8*
9#include "implicit.fh"
10#include "errquit.fh"
11#include "ciinfo.fh"
12      dimension hij(*),ind(4),indbar(*),iocc(*),
13     $     w1(nfmax,nfmax,nsmax),w2(nfmax2,nfmax,nsmax-1),
14     $     work1(*),work2(*),g(*),int12(*),int34(*)
15#include "stfunc.fh"
16c
17c     |I> is an orbital double replacement from |J>.
18c     Compute the hamiltonian matrix element in hij.
19c
20c     icase=1
21c        ii -> jj , ind(1)=i, ind(2)= j
22c        hij(u,v) = (ij|ij) delta(u,v) (u,v label spin functions)
23c     icase=2
24c        ik -> jj , ind(1)=i, ind(2)=k, ind(3)=j
25c        hij = (ij|kj)<u|Eij,kj|v> ... can be simplified
26c     icase=3
27c        ii -> jl , ind(1)=i, ind(2)=j, ind(3)=l
28c        hij = (ij|il)<u|Eij,il|v> ... can be simplified
29c     icase=4
30c        ik -> jl , ind(1)=i, ind(2)=k, ind(3)=j, ind(4)=l
31c        hij = (ij|kl)<u|Eij,kl|v> + (il|kj)<u|Eil,kj|v>
32c
33      numf = nf(ns)
34      goto (10,20,30,40) icase
35      call errquit('makhdb: invalid icase',icase, INPUT_ERR)
36c
37c     icase=1
38c        ii -> jj , ind(1)=i, ind(2)= j
39c        hij(u,v) = (ij|ij) delta(u,v) (u,v label spin functions)
40c
41 10   continue
42      numf2 = nf(ns)
43      call dfill(numf*numf,0.0d0,hij,1)
44      ij = itrian(ind(1),ind(2))
45      call dfill(numf,g(int12(ij)+int34(ij)),hij,numf+1)
46      return
47c
48c     icase=2
49c        ik -> jj , ind(1)=i, ind(2)=k, ind(3)=j
50c        hij = (ij|kj)<u|Eij,kj|v>
51c
52 20   continue
53      call selci_eijkj(hij,ind(1),ind(2),ind(3),ns,indbar,iocc,w1,w2,
54     $     work1,numf,numf2)
55      gg = g(intadr(itrian(ind(1),ind(3)),itrian(ind(2),ind(3))))
56      call dscal(numf*numf2,gg,hij,1)
57      return
58c
59c     icase=3
60c        ii -> jl , ind(1)=i, ind(2)=j, ind(3)=l
61c        hij = (ij|il)<u|Eij,il|v> ... can be simplified
62c
63 30   continue
64      call selci_eijil(hij,ind(1),ind(2),ind(3),ns,indbar,iocc,w1,w2,
65     $     work1,numf,numf2)
66      gg = g(intadr(itrian(ind(1),ind(2)),itrian(ind(1),ind(3))))
67      call dscal(numf*numf2,gg,hij,1)
68      return
69c
70c     icase=4
71c        ik -> jl , ind(1)=i, ind(2)=k, ind(3)=j, ind(4)=l
72c        hij = (ij|kl)<u|Eij,kl|v> + (il|kj)<u|Eil,kj|v>
73c
74 40   continue
75      call selci_eijkl(hij,ind(1),ind(3),ind(2),ind(4),ns,indbar,
76     $     iocc,w1,w2,work1,numf,numf2)
77      gg = g(intadr(itrian(ind(1),ind(3)),itrian(ind(2),ind(4))))
78      call dscal(numf*numf2,gg,hij,1)
79      call selci_eijkl(work2,ind(1),ind(4),ind(2),ind(3),ns,indbar,iocc,
80     $     w1,w2,work1,numf,numf2)
81      gg = g(intadr(itrian(ind(1),ind(4)),itrian(ind(3),ind(2))))
82      call daxpy(numf*numf2,gg,work2,1,hij,1)
83c
84      end
85C>
86C> @}
87