1      SUBROUTINE cr_eomccsd_tasim(d_x1,k_x1_offset,d_x2,k_x2_offset,
2     1                        d_t1,k_t1_offset,d_t2,k_t2_offset,
3     2                        d_f1,k_f1_offset,d_v2,k_v2_offset,
4     3                        d_e,k_e_offset,
5     4                        d_ex1,k_ex1_offset,size_ex1,
6     5                        d_ex2,k_ex2_offset,size_ex2,
7     6                        d_c1,k_c1_offset,size_c1,
8     7                        d_c2,k_c2_offset,size_c2,
9     8                        excit,energy1,energy2,size_t1,size_x1)
10C
11C     $Id$
12C
13c from this point on we assume that corresponding one- and two-
14c body components of R(d_x1,d_x2) and T(d_t1,d_t2) operators,
15c and corresponding excitation energy are available.
16c objects c1 and c2 are created in tce_energy right before
17c calling this procedure. The same structure of x1 and c1 and
18c x2 and c2 is assumed (irrep_c = irrep_x)
19c
20c 1. calculate R0
21c 2. calculate c1 and c2 vectors (c1excit,c2excit)
22c 3. calculate d
23c
24c local representations go only to "triple" procedures
25c see size_t1 size_x1  in the header of this procedure
26c
27      IMPLICIT NONE
28#include "global.fh"
29#include "mafdecls.fh"
30#include "util.fh"
31#include "errquit.fh"
32#include "tce.fh"
33#include "tce_main.fh"
34      integer d_t1
35      integer k_t1_offset
36      integer d_t2
37      integer k_t2_offset
38c --- sizes ---
39      integer size_c1,size_c2
40      integer size_ex1,size_ex2
41      integer size_d2
42c -------------
43      integer d_x1
44      integer k_x1_offset
45      integer d_x2
46      integer k_x2_offset
47      double precision r0xx  ! r0
48      double precision dr0xx ! r0*r0
49      double precision d1xx  ! <(singles)+singles>
50      double precision d2xx  ! <(doubles)+doubles>
51      double precision d12xx ! <(singles+doubles)+(singles+doubles)>
52      double precision d1xxt,d2xxt
53      double precision d1xxr,d2xxr
54      double precision d1xxtr,d2xxtr
55      double precision excit ! eomsd excitation energy
56      logical lr0            ! (true) r0*M3-calculated
57      integer d_ex1
58      integer k_ex1_offset
59      integer d_ex2
60      integer k_ex2_offset
61      integer d_c1
62      integer k_c1_offset
63      integer d_c2
64      integer k_c2_offset
65      integer d_d2
66      integer l_d2_offset
67      integer k_d2_offset
68c ----
69      integer d_f1
70      integer k_f1_offset
71      integer d_v2
72      integer k_v2_offset
73      integer d_e
74      integer k_e_offset
75      integer t_h1b, t_h1
76      integer t_h2b, t_h2
77      integer t_h3b, t_h3
78      integer t_p4b, t_p4
79      integer t_p5b, t_p5
80      integer t_p6b, t_p6
81      integer k_den,l_den
82      integer k_right1,l_right1  !r0*M3(T1,T2)
83      integer k_right2,l_right2  !M3(T1,T2,X1,X2)
84      integer k_right,l_right
85c --- cr_ccsd_t_E ---
86      integer k_left,l_left
87c -------------------
88c - T1/X1 LOCALIZATION -------------------
89      integer l_t1_local,k_t1_local
90      integer l_x1_local,k_x1_local
91      integer size_t1,size_x1
92c ---------------------------------------
93      integer size,i
94      integer g_energy
95      integer nxtask
96      integer next
97      integer nprocs
98      integer count
99c --- new intermediates ---
100c  cr_ccsd_t_N1 or cr_ccsd_t_N
101      integer d_i1_1,d_i1_2,d_i1_3
102      integer k_i1_offset_1,k_i1_offset_2,k_i1_offset_3
103      integer l_i1_offset_1,l_i1_offset_2,l_i1_offset_3
104c  cr_ccsd_t_E
105      integer d_i1_4,k_i1_offset_4,l_i1_offset_4
106c  cr_ccsd_t_N2
107      integer d_i2_1,d_i2_2,d_i2_3,d_i2_4,d_i2_5,d_i2_6
108      integer k_i2_offset_1,k_i2_offset_2,k_i2_offset_3
109      integer k_i2_offset_4,k_i2_offset_5,k_i2_offset_6
110      integer l_i2_offset_1,l_i2_offset_2,l_i2_offset_3
111      integer l_i2_offset_4,l_i2_offset_5,l_i2_offset_6
112c  q3rexpt
113      integer d_i3_1,d_i3_2
114      integer k_i3_offset_1,k_i3_offset_2
115      integer l_i3_offset_1,l_i3_offset_2
116c ----------------------------
117      logical is_active_t
118c ----------------------------
119      double precision energy1,energy2
120      double precision factor
121      double precision den1,num1
122      double precision den2,num2
123      double precision denex
124      character*255 filename
125      external nxtask
126c
127c Getting R0
128c
129      call  nr0_act(d_f1,d_e,d_t1,d_v2,d_x1,d_x2,k_f1_offset,
130     &k_e_offset,k_t1_offset,k_v2_offset,k_x1_offset,k_x2_offset)
131      call reconcilefile(d_e,1)
132      call get_block(d_e,r0xx,1,0)
133      if(dabs(excit).gt.1.0d-7) then
134         r0xx = r0xx/excit
135      else
136         write(6,1000)
137      end if
138      dr0xx = r0xx*r0xx
139      lr0 = .true.
140      if(dabs(r0xx).lt.1.0d-7) lr0 = .false.
141c
142c Now on ga with handle d_e we store corresponding R0 value
143c
144c
145c
146c Calculating one- and two-body overlaps
147c
148      if(lr0) then !symmetry of the reference
149         call tce_zero(d_ex1,size_ex1)
150         call tce_zero(d_ex2,size_ex2)
151         call c1_c1(d_t1,d_ex1,k_t1_offset,k_ex1_offset)
152         call reconcilefile(d_ex1,1)
153         d1xxt = 0.0d0
154         call get_block(d_ex1,d1xxt,1,0)
155         call t2t12_act(d_c2,d_t1,d_t2,k_c2_offset,
156     &        k_t1_offset,k_t2_offset)
157         call c2_c2_act(d_c2,d_ex2,k_c2_offset,k_ex2_offset)
158         call reconcilefile(d_ex2,1)
159         d2xxt = 0.0d0
160         call get_block(d_ex2,d2xxt,1,0)
161c forming vector R1T1+R2, d_d2's irrep corresponds to irrep_x
162c (in this case this is fully symmetric situation)
163c on d_c2 we have (T2+1/2T1*T1)|Phi\rangle
164            irrep_d=irrep_c
165            call tce_zero(d_ex1,size_ex1)
166            call tce_x2a_offset(l_d2_offset,k_d2_offset,size_d2)
167            call tce_filename('d2',filename)
168            call createfile(filename,d_d2,size_d2)
169            call tce_zero(d_d2,size_d2)
170         call c2excit2_act(d_d2,d_t1,d_x1,d_x2,k_d2_offset,
171     &     k_t1_offset,k_x1_offset,k_x2_offset)
172c <C2+ D2>  here
173      call c2_d2_act(d_c2,d_d2,d_ex1,k_c2_offset,k_d2_offset,
174     &     k_ex1_offset)
175         call reconcilefile(d_ex1,1)
176         d2xxtr = 0.0d0
177         call get_block(d_ex1,d2xxtr,1,0)
178            call deletefile(d_d2)
179            if (.not.ma_pop_stack(l_d2_offset))
180     1        call errquit("tce_energy: MA problem",36,MA_ERR)
181c
182         call tce_zero(d_ex1,size_ex1)
183         call tce_zero(d_ex2,size_ex2)
184         call x1_t1(d_ex1,d_t1,d_x1,k_ex1_offset,k_t1_offset,
185     &   k_x1_offset)
186         call reconcilefile(d_ex1,1)
187         d1xxtr=0.0d0
188         call get_block(d_ex1,d1xxtr,1,0)
189c
190         call tce_zero(d_ex1,size_ex1)
191         call tce_zero(d_ex2,size_ex2)
192         call tce_zero(d_c1,size_c1)
193         call tce_zero(d_c2,size_c2)
194         call c1excit2(d_c1,d_x1,k_c1_offset,k_x1_offset)
195         call c2excit2_act(d_c2,d_t1,d_x1,d_x2,k_c2_offset,
196     &     k_t1_offset,k_x1_offset,k_x2_offset)
197         call c1_c1(d_c1,d_ex1,k_c1_offset,k_ex1_offset)
198         call reconcilefile(d_ex1,1)
199         d1xxr = 0.0d0
200         call get_block(d_ex1,d1xxr,1,0)
201         call c2_c2_act(d_c2,d_ex2,k_c2_offset,k_ex2_offset)
202         call reconcilefile(d_ex2,1)
203         d2xxr = 0.0d0
204         call get_block(d_ex2,d2xxr,1,0)
205         d12xx = d1xxt*r0xx*r0xx+d2xxt*r0xx*r0xx+d1xxr+d2xxr+
206     &           2.0d0*d1xxtr*r0xx+2.0d0*d2xxtr*r0xx
207      else         !symmetry different form the symmetry of the reference
208c        <(R1+R0*T1)^{\dagger} (R1+R0*T1)>
209         call c1excit2(d_c1,d_x1,k_c1_offset,k_x1_offset)
210         call c2excit2_act(d_c2,d_t1,d_x1,d_x2,k_c2_offset,
211     &     k_t1_offset,k_x1_offset,k_x2_offset)
212         call c1_c1(d_c1,d_ex1,k_c1_offset,k_ex1_offset)
213         call reconcilefile(d_ex1,1)
214         d1xx = 0.0d0
215         call get_block(d_ex1,d1xx,1,0)
216c        <(R2+R1T1)^{\dagger} (R2+R1T1)>
217         call c2_c2_act(d_c2,d_ex2,k_c2_offset,k_ex2_offset)
218         call reconcilefile(d_ex2,1)
219         d2xx = 0.0d0
220         call get_block(d_ex2,d2xx,1,0)
221         d12xx = d1xx+d2xx
222      end if
223c - T1/X1 LOCALIZATION ----------
224c    opening l_t1_local and l_x1_local
225        if (.not.MA_PUSH_GET(mt_dbl,size_t1,'t1_local',
226     1      l_t1_local,k_t1_local))
227     1      call errquit('t1_local',1,MA_ERR)
228        if (.not.MA_PUSH_GET(mt_dbl,size_x1,'x1_local',
229     1      l_x1_local,k_x1_local))
230     1      call errquit('x1_local',1,MA_ERR)
231        call ma_zero(dbl_mb(k_t1_local),size_t1)
232        call ma_zero(dbl_mb(k_x1_local),size_x1)
233c    copy d_t1 ==> l_t1_local
234c    copy x1(ivec) ==> l_x1_local
235        call ga_get(d_t1,1,size_t1,1,1,dbl_mb(k_t1_local),1)
236        call ga_get(d_x1,1,size_x1,1,1,dbl_mb(k_x1_local),1)
237c -------------------------------
238c
239c     Caution! k_right & k_den are not even allocated yet
240c     but they won't be used.
241c --------------- initialization -----------------------------------
242c
243cc      if(lr0) then
244cc      call cr_ccsd_t_N_act(dbl_mb(k_right),d_f1,d_i1_1,d_i1_2,
245cc     1  k_t1_local,d_t2,d_v2,k_f1_offset,k_i1_offset_1,k_i1_offset_2,
246cc     2  k_t1_offset,k_t2_offset,k_v2_offset,l_i1_offset_1,
247cc     3  l_i1_offset_2,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,1)
248cc     end if
249c
250cc      call creomsd_tact_mem(dbl_mb(k_right),d_f1,d_i2_1,d_i2_2,
251cc     &d_i2_3,d_i2_4,k_t1_local,d_t2,d_v2,k_x1_local,
252cc     &d_x2,k_f1_offset,k_i2_offset_1,
253cc     &k_i2_offset_2,k_i2_offset_3,k_i2_offset_4,k_t1_offset,k_t2_offset,
254cc     &k_v2_offset,k_x1_offset,k_x2_offset,l_i2_offset_1,l_i2_offset_2,
255cc     &l_i2_offset_3,l_i2_offset_4,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,1)
256c
257c
258      call cr_eomt_simac(dbl_mb(k_right),d_i2_1,d_i2_2,d_t2,d_v2,
259     &d_x1,d_x2,k_i2_offset_1,k_i2_offset_2,k_t2_offset,k_v2_offset,
260     &k_x1_offset,k_x2_offset,l_i2_offset_1,l_i2_offset_2,
261     &t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,1)
262c
263c
264      if(lr0) then
265      call cr_ccsd_t_E_act(dbl_mb(k_left),d_i1_4,
266     1  k_t1_local,d_t2,k_i1_offset_4,k_t1_offset,k_t2_offset,
267     2  l_i1_offset_4,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,1)
268      call q3rexpt2_act(dbl_mb(k_left),d_i3_1,
269     &k_t1_local,d_t2,k_x1_local,d_x2,
270     &k_i3_offset_1,k_t1_offset,k_t2_offset,k_x1_offset,k_x2_offset,
271     &l_i3_offset_1,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,1)
272      else
273      call q3rexpt2_act(dbl_mb(k_left),d_i3_1,k_t1_local,d_t2,
274     &k_x1_local,d_x2,
275     &k_i3_offset_1,k_t1_offset,k_t2_offset,k_x1_offset,k_x2_offset,
276     &l_i3_offset_1,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,1)
277      end if
278c ------------------------------------------------------------------
279c
280c     Get the numerator
281c
282      num1 = 0.0d0
283      den1 = 0.0d0
284c
285      num2 = 0.0d0
286      den2 = 0.0d0
287c
288      if (.not.ga_create(mt_dbl,1,1,'perturbative',1,1,g_energy))
289     1  call errquit('ccsd_t: GA problem',0,GA_ERR)
290      nprocs = GA_NNODES()
291      count = 0
292      next = nxtask(nprocs,1)
293      do t_p4b = noab+1,noab+nvab
294       do t_p5b = t_p4b,noab+nvab
295        do t_p6b = t_p5b,noab+nvab
296         do t_h1b = 1,noab
297          do t_h2b = t_h1b,noab
298           do t_h3b = t_h2b,noab
299            if (next.eq.count) then
300            IF(is_active_t(t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b)) THEN
301            if (int_mb(k_spin+t_p4b-1)
302     1         +int_mb(k_spin+t_p5b-1)
303     2         +int_mb(k_spin+t_p6b-1)
304     3      .eq.int_mb(k_spin+t_h1b-1)
305     4         +int_mb(k_spin+t_h2b-1)
306     5         +int_mb(k_spin+t_h3b-1)) then
307            if ((.not.restricted).or.
308     1         (int_mb(k_spin+t_p4b-1)
309     1         +int_mb(k_spin+t_p5b-1)
310     2         +int_mb(k_spin+t_p6b-1)
311     3         +int_mb(k_spin+t_h1b-1)
312     4         +int_mb(k_spin+t_h2b-1)
313     5         +int_mb(k_spin+t_h3b-1).le.8)) then
314            if (ieor(int_mb(k_sym+t_p4b-1),
315     1          ieor(int_mb(k_sym+t_p5b-1),
316     2          ieor(int_mb(k_sym+t_p6b-1),
317     3          ieor(int_mb(k_sym+t_h1b-1),
318     4          ieor(int_mb(k_sym+t_h2b-1),
319     5               int_mb(k_sym+t_h3b-1)))))).eq.irrep_x) then
320c
321c Symmetry control (above)
322c
323            size = int_mb(k_range+t_p4b-1)
324     1           * int_mb(k_range+t_p5b-1)
325     2           * int_mb(k_range+t_p6b-1)
326     3           * int_mb(k_range+t_h1b-1)
327     4           * int_mb(k_range+t_h2b-1)
328     5           * int_mb(k_range+t_h3b-1)
329            if (.not.MA_PUSH_GET(mt_dbl,size,'moment 2,3',
330     1        l_right,k_right)) call errquit('eomccsd_t',3,MA_ERR)
331            if (.not.MA_PUSH_GET(mt_dbl,size,'denominator',
332     1        l_left,k_left)) call errquit('ccsd_t',3,MA_ERR)
333ccx            do i = 1, size
334ccx             dbl_mb(k_right+i-1) = 0.0d0
335ccx             dbl_mb(k_left+i-1) = 0.0d0
336ccx            enddo
337c zeroing ---
338        call dfill(size, 0.0d0, dbl_mb(k_right), 1)
339        call dfill(size, 0.0d0, dbl_mb(k_left), 1)
340c -----------
341c
342c Moments are calculated here
343c
344c
345      if(lr0) then
346c
347cc      call cr_ccsd_t_N_act(dbl_mb(k_right),d_f1,d_i1_1,d_i1_2,
348cc     1  k_t1_local,d_t2,d_v2,k_f1_offset,k_i1_offset_1,k_i1_offset_2,
349cc     2  k_t1_offset,k_t2_offset,k_v2_offset,l_i1_offset_1,
350cc     3  l_i1_offset_2,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,2)
351c
352      call ccsd_t_d_act(dbl_mb(k_right),d_t2,d_v2,k_t2_offset,
353     &k_v2_offset,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,2)
354c
355cc            do i = 1, size
356cc             dbl_mb(k_right+i-1) = r0xx*dbl_mb(k_right+i-1)
357cc            enddo
358        call dscal(size,r0xx,dbl_mb(k_right),1)
359c
360      end if
361c
362cc      call creomsd_tact_mem(dbl_mb(k_right),d_f1,d_i2_1,d_i2_2,
363cc     &d_i2_3,d_i2_4,k_t1_local,d_t2,d_v2,
364cc     &k_x1_local,d_x2,k_f1_offset,k_i2_offset_1,
365cc     &k_i2_offset_2,k_i2_offset_3,k_i2_offset_4,k_t1_offset,k_t2_offset,
366cc     &k_v2_offset,k_x1_offset,k_x2_offset,l_i2_offset_1,l_i2_offset_2,
367cc     &l_i2_offset_3,l_i2_offset_4,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,2)
368c
369      call cr_eomt_simac(dbl_mb(k_right),d_i2_1,d_i2_2,d_t2,d_v2,
370     &d_x1,d_x2,k_i2_offset_1,k_i2_offset_2,k_t2_offset,k_v2_offset,
371     &k_x1_offset,k_x2_offset,l_i2_offset_1,l_i2_offset_2,
372     &t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,2)
373c
374c Q3(R0+R1+R2)exp(T1+T2)|Ref> calculated here
375c
376      if(lr0) then
377      call cr_ccsd_t_E_act(dbl_mb(k_left),d_i1_4,
378     1  k_t1_local,d_t2,k_i1_offset_4,k_t1_offset,k_t2_offset,
379     2  l_i1_offset_4,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,2)
380c
381cc            do i = 1, size
382cc             dbl_mb(k_left+i-1) = r0xx*dbl_mb(k_left+i-1)
383cc            enddo
384        call dscal(size,r0xx,dbl_mb(k_left),1)
385c
386      call q3rexpt2_act(dbl_mb(k_left),d_i3_1,k_t1_local,d_t2,
387     &k_x1_local,d_x2,
388     &k_i3_offset_1,k_t1_offset,k_t2_offset,k_x1_offset,k_x2_offset,
389     &l_i3_offset_1,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,2)
390      else
391      call q3rexpt2_act(dbl_mb(k_left),d_i3_1,k_t1_local,d_t2,
392     &k_x1_local,d_x2,
393     &k_i3_offset_1,k_t1_offset,k_t2_offset,k_x1_offset,k_x2_offset,
394     &l_i3_offset_1,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,2)
395      end if
396c
397            if (restricted) then
398              factor = 2.0d0
399            else
400              factor = 1.0d0
401            endif
402            if ((t_p4b.eq.t_p5b).and.(t_p5b.eq.t_p6b)) then
403              factor = factor / 6.0d0
404            else if ((t_p4b.eq.t_p5b).or.(t_p5b.eq.t_p6b)) then
405              factor = factor / 2.0d0
406            endif
407            if ((t_h1b.eq.t_h2b).and.(t_h2b.eq.t_h3b)) then
408              factor = factor / 6.0d0
409            else if ((t_h1b.eq.t_h2b).or.(t_h2b.eq.t_h3b)) then
410              factor = factor / 2.0d0
411            endif
412c
413c
414c
415            i = 0
416            do t_p4 = 1, int_mb(k_range+t_p4b-1)
417             do t_p5 = 1, int_mb(k_range+t_p5b-1)
418              do t_p6 = 1, int_mb(k_range+t_p6b-1)
419               do t_h1 = 1, int_mb(k_range+t_h1b-1)
420                do t_h2 = 1, int_mb(k_range+t_h2b-1)
421                 do t_h3 = 1, int_mb(k_range+t_h3b-1)
422                  i = i + 1
423            denex=-dbl_mb(k_evl_sorted+int_mb(k_offset+t_p4b-1)+t_p4-1)
424     4         -dbl_mb(k_evl_sorted+int_mb(k_offset+t_p5b-1)+t_p5-1)
425     5         -dbl_mb(k_evl_sorted+int_mb(k_offset+t_p6b-1)+t_p6-1)
426     6         +dbl_mb(k_evl_sorted+int_mb(k_offset+t_h1b-1)+t_h1-1)
427     7         +dbl_mb(k_evl_sorted+int_mb(k_offset+t_h2b-1)+t_h2-1)
428     8         +dbl_mb(k_evl_sorted+int_mb(k_offset+t_h3b-1)+t_h3-1)
429     9         +excit
430c numerator
431                  num1 = num1 + factor *
432     1            (dbl_mb(k_right+i-1)*dbl_mb(k_right+i-1))/denex
433                  num1 = num1 + factor *
434     1            dbl_mb(k_left+i-1)*dbl_mb(k_right+i-1)
435c denominator
436c
437                  den1 = den1 + factor *
438     1            (dbl_mb(k_left+i-1)*dbl_mb(k_right+i-1))/denex
439                  den1 = den1 + factor *
440     1            (dbl_mb(k_left+i-1)*dbl_mb(k_left+i-1))
441c
442                 enddo
443                enddo
444               enddo
445              enddo
446             enddo
447            enddo
448c
449c
450c
451            if (.not.MA_POP_STACK(l_left))
452     1        call errquit('eomccsd_t',6,MA_ERR)
453            if (.not.MA_POP_STACK(l_right))
454     1        call errquit('eomccsd_t',6,MA_ERR)
455            endif
456            endif
457            endif
458            endif
459            next = nxtask(nprocs,1)
460            endif
461            count = count + 1
462           enddo
463          enddo
464         enddo
465        enddo
466       enddo
467      enddo
468      next = nxtask(-nprocs,1)
469c --- toggle = 3 ---
470      if(lr0) then
471      call q3rexpt2_act(dbl_mb(k_left),d_i3_1,k_t1_local,d_t2,
472     &k_x1_local,d_x2,
473     &k_i3_offset_1,k_t1_offset,k_t2_offset,k_x1_offset,k_x2_offset,
474     &l_i3_offset_1,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,3)
475      call cr_ccsd_t_E_act(dbl_mb(k_left),d_i1_4,
476     1  k_t1_local,d_t2,k_i1_offset_4,k_t1_offset,k_t2_offset,
477     2  l_i1_offset_4,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,3)
478      else
479      call q3rexpt2_act(dbl_mb(k_left),d_i3_1,k_t1_local,d_t2,
480     &k_x1_local,d_x2,
481     &k_i3_offset_1,k_t1_offset,k_t2_offset,k_x1_offset,k_x2_offset,
482     &l_i3_offset_1,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,3)
483      end if
484c
485cc      call creomsd_tact_mem(dbl_mb(k_right),d_f1,d_i2_1,d_i2_2,
486cc     &d_i2_3,d_i2_4,k_t1_local,d_t2,d_v2,
487cc     &k_x1_local,d_x2,k_f1_offset,k_i2_offset_1,
488cc     &k_i2_offset_2,k_i2_offset_3,k_i2_offset_4,k_t1_offset,k_t2_offset,
489cc     &k_v2_offset,k_x1_offset,k_x2_offset,l_i2_offset_1,l_i2_offset_2,
490cc     &l_i2_offset_3,l_i2_offset_4,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,3)
491c
492      call cr_eomt_simac(dbl_mb(k_right),d_i2_1,d_i2_2,d_t2,d_v2,
493     &d_x1,d_x2,k_i2_offset_1,k_i2_offset_2,k_t2_offset,k_v2_offset,
494     &k_x1_offset,k_x2_offset,l_i2_offset_1,l_i2_offset_2,
495     &t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,3)
496c
497cc      if(lr0) then
498cc      call cr_ccsd_t_N_act(dbl_mb(k_right),d_f1,d_i1_1,d_i1_2,
499cc     1  k_t1_local,d_t2,d_v2,k_f1_offset,k_i1_offset_1,k_i1_offset_2,
500cc     2  k_t1_offset,k_t2_offset,k_v2_offset,l_i1_offset_1,
501cc     3  l_i1_offset_2,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,3)
502cc      end if
503c
504c ------------------
505c
506      call ga_zero(g_energy)
507      call ga_acc(g_energy,1,1,1,1,num1,1,1.0d0)
508      call ga_sync()
509      call ga_get(g_energy,1,1,1,1,num1,1)
510c
511      call ga_zero(g_energy)
512      call ga_acc(g_energy,1,1,1,1,den1,1,1.0d0)
513      call ga_sync()
514      call ga_get(g_energy,1,1,1,1,den1,1)
515c
516      call ga_zero(g_energy)
517      call ga_acc(g_energy,1,1,1,1,num2,1,1.0d0)
518      call ga_sync()
519      call ga_get(g_energy,1,1,1,1,num2,1)
520c
521      call ga_zero(g_energy)
522      call ga_acc(g_energy,1,1,1,1,den2,1,1.0d0)
523      call ga_sync()
524      call ga_get(g_energy,1,1,1,1,den2,1)
525c
526      if (.not.ga_destroy(g_energy))
527     1  call errquit('creom_t: GA problem',1,GA_ERR)
528       energy1 = num1/(dr0xx+d12xx+den1)
529      call util_flush(6)
530c - T1/X1 LOCALIZATION ------
531         if(.not.MA_POP_STACK(l_x1_local))
532     &      call errquit('l_x1_local',4,MA_ERR)
533         if(.not.MA_POP_STACK(l_t1_local))
534     &      call errquit('l_t1_local',4,MA_ERR)
535c ---------------------------
536 1000    format('corresponding excitation energy = 0')
537c
538      return
539      end
540c
541
542