1      subroutine rohf_hessv_1e3( basis, geom, nmo, nclosed, nopen,
2     $     g_fcv, g_fpv, g_fcp,
3     $     g_x, g_ax, omega)
4C     $Id$
5c
6c ... jochen: modified version of rohf_hessv_1e that goes along with
7c     rohf_hessv_2e3 and the other routines used by cphf_solve3
8c
9      implicit none
10#include "errquit.fh"
11#include "mafdecls.fh"
12#include "global.fh"
13c
14      integer basis, geom
15      integer nmo, nclosed, nopen
16      integer g_fcv, g_fpv, g_fcp
17      integer g_x
18      integer g_ax
19      double precision omega
20c
21      integer nvir, voff, ooff, oend, xoff, xend
22      integer g_tmp1, g_tmp2, g_tmp3
23      double precision four, mfour, one, mone, two, mtwo, zero
24      parameter (four=4.0d0, mfour=-4.0d0, one=1.0d0, mone=-1.0d0,
25     $   two=2.0d0, mtwo=-2.0d0, zero=0.0d0)
26c
27      integer ivec, nvec, gtype, vlen
28c
29      call ga_inquire(g_x,gtype,vlen,nvec)
30c
31      nvir = nmo - nclosed - nopen
32      voff = nclosed + nopen + 1
33      ooff = nclosed + 1
34      oend = nclosed + nopen
35c
36c     Reshape argument vector by copying patches into matrix
37c
38c      cv      cv  cv       cv  cv
39c     S   =  4.F  .x   -  4.x  .F
40c
41      if (.not. ga_create(MT_DBL, nmo, nmo, 'rohf_hv1e: tmp1',
42     $     32, 32, g_tmp1)) call errquit('rohf_hv1e: tmp1', 0, GA_ERR)
43      if (.not. ga_create(MT_DBL, nmo, nmo, 'rohf_hv1e: tmp2',
44     $     32, 32, g_tmp2)) call errquit('rohf_hv1e: tmp2', 0, GA_ERR)
45      if (nopen .gt. 0) then
46         if (.not. ga_create(MT_DBL, nmo, nmo, 'rohf_hv1e: tmp3',
47     $        32, 32, g_tmp3)) call errquit('rohf_hv1e: tmp3', 0,
48     &       GA_ERR)
49      endif
50c
51      do ivec = 1, nvec
52        call ga_zero(g_tmp1)
53c
54c       copy input vector g_x in temp array g_tmp1
55c       (g_x holds the elements of the A-matrix)
56c
57        call ga_vec_to_mat(g_tmp1, voff, nmo, 1, nclosed,
58     $     g_x(ipm), 1, ivec)
59*       call ga_copy_patch('n', g_x, 1, (nvir*nclosed), ivec, ivec,
60*       $        g_tmp1, voff, nmo, 1, nclosed )
61        call ga_zero(g_tmp2)
62        call ga_matmul_patch('n','n',four,zero,
63     $     g_fcv, voff, nmo, voff, nmo, ! virtual part of Fock-operator
64     $     g_tmp1, voff, nmo, 1, nclosed, ! input vector, virt-occ
65     $     g_tmp2, voff, nmo, 1, nclosed ) ! 4 times product
66        call ga_matmul_patch('n','n',mfour,one,
67     $     g_tmp1, voff, nmo, 1, nclosed,  ! input vector, virt-occ
68     $     g_fcv, 1, nclosed, 1, nclosed,  ! occ. part of Fock-operator
69     $     g_tmp2, voff, nmo, 1, nclosed ) ! -4 times add to prev. res.
70c
71c ... jochen: if this makes any sense so far, we will now
72c       add 4 times omega to this.
73c
74       call ga_add_patch(
75     &     1d0, g_tmp2, voff, nmo, 1, nclosed,
76     &     four * omega, g_tmp1, voff, nmo, 1, nclosed,
77     &     g_tmp2, voff, nmo, 1, nclosed,
78c
79        if (nopen .gt. 0) then
80c
81          xoff = nclosed*nvir + 1
82          xend = (nclosed + nopen)*nvir
83          call ga_copy_patch('n', g_x, xoff, xend, ivec, ivec,
84     $       g_tmp1, voff, nmo, ooff, oend )
85          xoff = nvir*(nclosed+nopen) + 1
86          xend = nvir*(nclosed+nopen) + nclosed*nopen
87          call ga_copy_patch('n', g_x, xoff, xend, ivec, ivec,
88     $       g_tmp1, ooff, oend, 1, nclosed )
89c
90c       pv       pv  pv       pv  pv
91c     S   =  2.F  .x   -  2.x  .F
92c
93            call ga_matmul_patch('n','n',two,one,
94     $           g_fpv, voff, nmo, voff, nmo,
95     $           g_tmp1, voff, nmo, ooff, oend,
96     $           g_tmp2, voff, nmo, ooff, oend )
97            call ga_matmul_patch('n','n', mtwo, one,
98     $           g_tmp1, voff, nmo, ooff, oend,
99     $           g_fpv, ooff, oend, ooff, oend,
100     $           g_tmp2, voff, nmo, ooff, oend )
101c
102c       cp       cp  cp       cp  cp
103c     S   =  2.F  .x   -  2.x  .F
104c
105            call ga_matmul_patch('n','n',two,one,
106     $           g_fcp, ooff, oend, ooff, oend,
107     $           g_tmp1, ooff, oend, 1, nclosed,
108     $           g_tmp2, ooff, oend, 1, nclosed )
109            call ga_matmul_patch('n','n',mtwo,one,
110     $           g_tmp1, ooff, oend, 1, nclosed,
111     $           g_fcp, 1, nclosed, 1, nclosed,
112     $           g_tmp2, ooff, oend, 1, nclosed )
113c
114c       cv         cp     pv   cp
115c     S    += (2.F   +  F  ).x
116c
117c       cp         cp     pv   cv
118c     S    += (2.F   +  F  ).x
119c
120            call ga_dadd(two,g_fcp,one,g_fpv,g_tmp3)
121            call ga_matmul_patch('n', 'n', one, one,
122     $           g_tmp3, voff, nmo, ooff, oend,
123     $           g_tmp1, ooff, oend, 1, nclosed,
124     $           g_tmp2, voff, nmo, 1, nclosed )
125            call ga_matmul_patch('t', 'n', one, one,
126     $           g_tmp3, ooff, oend, voff, nmo,
127     $           g_tmp1, voff, nmo, 1, nclosed,
128     $           g_tmp2, ooff, oend, 1, nclosed )
129c
130c       cv         pv     cp   pv
131c     S    -= (2.F   +  F  ).x
132c
133c       pv         pv     cp   cv
134c     S    -= (2.F   +  F  ).x
135c
136            call ga_dadd(one,g_fcp,two,g_fpv,g_tmp3)
137            call ga_matmul_patch('n','n',mone,one,
138     $           g_tmp1, voff, nmo, ooff, oend,
139     $           g_tmp3, ooff, oend, 1, nclosed,
140     $           g_tmp2, voff, nmo, 1, nclosed )
141            call ga_matmul_patch('n','t',mone,one,
142     $           g_tmp1, voff, nmo, 1, nclosed,
143     $           g_tmp3, 1, nclosed, ooff, oend,
144     $           g_tmp2, voff, nmo, ooff, oend )
145c
146c       pv       cp     pv   cp
147c     S    += (F   -  F  ).x
148c
149c       cp       cp     pv   pv
150c     S    += (F   -  F  ).x
151c
152            call ga_dadd(one,g_fcp,mone,g_fpv,g_tmp3)
153            call ga_matmul_patch('n', 't', one, one,
154     $           g_tmp3, voff, nmo, 1, nclosed,
155     $           g_tmp1, 1, nclosed, ooff, oend,
156     $           g_tmp2, voff, nmo, ooff, oend )
157            call ga_matmul_patch('t', 'n', one, one,
158     $           g_tmp1, ooff, oend, voff, nmo,
159     $           g_tmp3, voff, nmo, 1, nclosed,
160     $           g_tmp2, ooff, oend, 1, nclosed )
161c
162         endif ! nopen.gt.0
163c
164         xend = nclosed*nvir
165c
166         call ga_mat_to_vec(g_tmp2, voff, nmo, 1, nclosed,
167     $       g_ax, 1, ivec, 1.0d0, '+')
168*         call ga_dadd_patch( one, g_tmp2, voff, nmo, 1, nclosed,
169*     $        one, g_ax, 1, xend, ivec, ivec,
170*     $        g_ax, 1, xend, ivec, ivec)
171c
172         if (nopen .gt. 0) then
173            xoff = xend + 1
174            xend = xend + nopen*nvir
175*
176            call ga_mat_to_vec(g_tmp2, voff, nmo, ooff, oend,
177     $         g_ax, xoff, ivec, 1.0d0, '+')
178*            call ga_dadd_patch( one, g_tmp2, voff, nmo, ooff, oend,
179*     $           one, g_ax, xoff, xend, ivec, ivec,
180*     $           g_ax, xoff, xend, ivec, ivec)
181            xoff = xend + 1
182            xend = xend + nopen*nclosed
183            call ga_mat_to_vec(g_tmp2, ooff, oend, 1, nclosed,
184     $         g_ax, xoff, ivec, 1.0d0, '+')
185*            call ga_dadd_patch( one, g_tmp2, ooff, oend, 1, nclosed,
186*     $           one, g_ax, xoff, xend, ivec, ivec,
187*     $           g_ax, xoff, xend, ivec, ivec)
188         endif
189      enddo
190c
191      if (.not. ga_destroy(g_tmp1)) call errquit('rohf_hv1?',0, GA_ERR)
192      if (.not. ga_destroy(g_tmp2)) call errquit('rohf_hv1?',0, GA_ERR)
193      if (nopen .gt. 0) then
194         if (.not. ga_destroy(g_tmp3)) call errquit('rohf_hv1?',0,
195     &       GA_ERR)
196      endif
197c
198      end
199