~ubuntu-branches/ubuntu/utopic/nwchem/utopic

« back to all changes in this revision

Viewing changes to src/tce/mrcc/ccsd_t/tce_mrcc_ccsdpt.F

  • Committer: Package Import Robot
  • Author(s): Michael Banck, Daniel Leidert, Andreas Tille, Michael Banck
  • Date: 2013-07-04 12:14:55 UTC
  • mfrom: (1.1.2)
  • Revision ID: package-import@ubuntu.com-20130704121455-5tvsx2qabor3nrui
Tags: 6.3-1
* New upstream release.
* Fixes anisotropic properties (Closes: #696361).
* New features include:
  + Multi-reference coupled cluster (MRCC) approaches
  + Hybrid DFT calculations with short-range HF 
  + New density-functionals including Minnesota (M08, M11) and HSE hybrid
    functionals
  + X-ray absorption spectroscopy (XAS) with TDDFT
  + Analytical gradients for the COSMO solvation model
  + Transition densities from TDDFT 
  + DFT+U and Electron-Transfer (ET) methods for plane wave calculations
  + Exploitation of space group symmetry in plane wave geometry optimizations
  + Local density of states (LDOS) collective variable added to Metadynamics
  + Various new XC functionals added for plane wave calculations, including
    hybrid and range-corrected ones
  + Electric field gradients with relativistic corrections 
  + Nudged Elastic Band optimization method
  + Updated basis sets and ECPs 

[ Daniel Leidert ]
* debian/watch: Fixed.

[ Andreas Tille ]
* debian/upstream: References

[ Michael Banck ]
* debian/upstream (Name): New field.
* debian/patches/02_makefile_flags.patch: Refreshed.
* debian/patches/06_statfs_kfreebsd.patch: Likewise.
* debian/patches/07_ga_target_force_linux.patch: Likewise.
* debian/patches/05_avoid_inline_assembler.patch: Removed, no longer needed.
* debian/patches/09_backported_6.1.1_fixes.patch: Likewise.
* debian/control (Build-Depends): Added gfortran-4.7 and gcc-4.7.
* debian/patches/10_force_gcc-4.7.patch: New patch, explicitly sets
  gfortran-4.7 and gcc-4.7, fixes test suite hang with gcc-4.8 (Closes:
  #701328, #713262).
* debian/testsuite: Added tests for COSMO analytical gradients and MRCC.
* debian/rules (MRCC_METHODS): New variable, required to enable MRCC methods.

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
ckbn ccsd_t.F modified to tce_mrcc_ccsdpt_uncoup_pt3.F 
 
2
 
 
3
       SUBROUTINE tce_mrcc_ccsdpt(d_t1,k_t1_offset,
 
4
     +                           d_t2,k_t2_offset,
 
5
     +                           d_v2,k_v2_offset,
 
6
     +                           d_f1,
 
7
     +                           k_f1_offset,
 
8
ckbn +                           energy1,energy2,energy3,
 
9
     +                           energy3,
 
10
     +                           iref)
 
11
c     1                           iref,nref)
 
12
 
 
13
      IMPLICIT NONE
 
14
#include "mafdecls.fh"
 
15
#include "tcgmsg.fh"
 
16
#include "global.fh"
 
17
#include "bas.fh"
 
18
#include "geom.fh"
 
19
#include "rtdb.fh"
 
20
#include "sym.fh"
 
21
#include "util.fh"
 
22
#include "msgids.fh"
 
23
#include "stdio.fh"
 
24
#include "sf.fh"
 
25
#include "inp.fh"
 
26
#include "errquit.fh"
 
27
#include "tce.fh"
 
28
#include "tce_main.fh"
 
29
#include "tce_hetio.fh"
 
30
#include "tce_diis.fh"
 
31
#include "tce_prop.fh"
 
32
#include "tce_restart.fh"
 
33
#include "tce_mrcc.fh"
 
34
 
 
35
c      integer nref
 
36
      integer iref
 
37
      integer d_t1
 
38
      integer k_t1_offset
 
39
      integer d_t2
 
40
      integer k_t2_offset
 
41
      integer d_v2
 
42
      integer k_v2_offset
 
43
      integer t_h1b, t_h1
 
44
      integer t_h2b, t_h2
 
45
      integer t_h3b, t_h3
 
46
      integer t_p4b, t_p4
 
47
      integer t_p5b, t_p5
 
48
      integer t_p6b, t_p6
 
49
      integer k_singles,l_singles
 
50
      integer k_doubles,l_doubles
 
51
ckbn-2
 
52
c      integer k_q3fnt2,l_q3fnt2
 
53
      integer size,i,j
 
54
      integer nxtask
 
55
      integer next
 
56
      integer nprocs
 
57
      integer count
 
58
      integer offset_p4,offset_p5,offset_p6
 
59
      integer offset_h1,offset_h2,offset_h3
 
60
      integer range_p4,range_p5,range_p6
 
61
      integer range_h1,range_h2,range_h3
 
62
      double precision energy(3)
 
63
ckbn      double precision energy1,energy2
 
64
ckbn-2
 
65
      double precision energy3
 
66
      double precision factor,denom
 
67
      double precision denom_p4,denom_p5,denom_p6
 
68
      double precision denom_h1,denom_h2,denom_h3
 
69
      external nxtask
 
70
ckbn-3
 
71
      logical nodezero
 
72
      integer d_f1,k_f1_offset
 
73
      double precision dnorm1,dnorm2,dmaxt1,dmaxt2
 
74
      double precision dmint1,dmint2
 
75
      integer ioff
 
76
      integer orbspin(6),orbindex(6),off(6),blck(6)
 
77
      INTEGER NXTASKsub
 
78
      EXTERNAL NXTASKsub
 
79
 
 
80
ckbn-2
 
81
      nodezero=(ga_nodeid().eq.0)
 
82
 
 
83
ckbn-2
 
84
c      if(nodezero) write(*,*)"Entering tce_mrcc_ccsdpt"
 
85
 
 
86
      nprocs = GA_NNODES()
 
87
      count = 0
 
88
      next = nxtask(nprocs,1)
 
89
 
 
90
 
 
91
 
 
92
      energy(1)=0.0d0
 
93
      energy(2)=0.0d0
 
94
      energy(3)=0.0d0
 
95
c      energy1  =0.0d0
 
96
c      energy2  =0.0d0
 
97
      energy3  =0.0d0
 
98
      do t_p4b = noab+1,noab+nvab
 
99
       range_p4 = int_mb(k_range+t_p4b-1)
 
100
       offset_p4 = k_evl_sorted+int_mb(k_offset+t_p4b-1)-1
 
101
       do t_p5b = t_p4b,noab+nvab
 
102
        range_p5 = int_mb(k_range+t_p5b-1)
 
103
        offset_p5 = k_evl_sorted+int_mb(k_offset+t_p5b-1)-1
 
104
        do t_p6b = t_p5b,noab+nvab
 
105
         range_p6 = int_mb(k_range+t_p6b-1)
 
106
         offset_p6 = k_evl_sorted+int_mb(k_offset+t_p6b-1)-1
 
107
         do t_h1b = 1,noab
 
108
          range_h1 = int_mb(k_range+t_h1b-1)
 
109
          offset_h1 = k_evl_sorted+int_mb(k_offset+t_h1b-1)-1
 
110
          do t_h2b = t_h1b,noab
 
111
           range_h2 = int_mb(k_range+t_h2b-1)
 
112
           offset_h2 = k_evl_sorted+int_mb(k_offset+t_h2b-1)-1
 
113
           do t_h3b = t_h2b,noab
 
114
            range_h3 = int_mb(k_range+t_h3b-1)
 
115
            offset_h3 = k_evl_sorted+int_mb(k_offset+t_h3b-1)-1
 
116
            if (next.eq.count) then                        
 
117
            if (int_mb(k_spin+t_p4b-1)
 
118
     1         +int_mb(k_spin+t_p5b-1)
 
119
     2         +int_mb(k_spin+t_p6b-1)
 
120
     3      .eq.int_mb(k_spin+t_h1b-1)
 
121
     4         +int_mb(k_spin+t_h2b-1)
 
122
     5         +int_mb(k_spin+t_h3b-1)) then
 
123
            if ((.not.restricted).or.
 
124
     1         (int_mb(k_spin+t_p4b-1)
 
125
     1         +int_mb(k_spin+t_p5b-1)
 
126
     2         +int_mb(k_spin+t_p6b-1)
 
127
     3         +int_mb(k_spin+t_h1b-1)
 
128
     4         +int_mb(k_spin+t_h2b-1)
 
129
     5         +int_mb(k_spin+t_h3b-1).le.8)) then
 
130
            if (ieor(int_mb(k_sym+t_p4b-1),
 
131
     1          ieor(int_mb(k_sym+t_p5b-1),
 
132
     2          ieor(int_mb(k_sym+t_p6b-1),
 
133
     3          ieor(int_mb(k_sym+t_h1b-1),
 
134
     4          ieor(int_mb(k_sym+t_h2b-1),
 
135
     5               int_mb(k_sym+t_h3b-1)))))).eq.0) then
 
136
 
 
137
      IF (.not.((log_mb(k_active_tmpm(iref)+t_p4b-1).eqv..true.).and.
 
138
     1  (log_mb(k_active_tmpm(iref)+t_p5b-1).eqv..true.).and.
 
139
     1  (log_mb(k_active_tmpm(iref)+t_p6b-1).eqv..true.).and.
 
140
     1  (log_mb(k_active_tmpm(iref)+t_h1b-1).eqv..true.).and.
 
141
     1  (log_mb(k_active_tmpm(iref)+t_h2b-1).eqv..true.).and.
 
142
     1  (log_mb(k_active_tmpm(iref)+t_h3b-1).eqv..true.))) then
 
143
 
 
144
            size = range_p4 * range_p5 * range_p6
 
145
     3           * range_h1 * range_h2 * range_h3
 
146
            if (.not.MA_PUSH_GET(mt_dbl,size,'(T) singles',l_singles,
 
147
     1        k_singles)) call errquit('ccsd_t: MA error',1,MA_ERR)
 
148
            if (.not.MA_PUSH_GET(mt_dbl,size,'(T) doubles',l_doubles,
 
149
     1        k_doubles)) call errquit('ccsd_t: MA error',2,MA_ERR)
 
150
            do i = 0, size-1
 
151
              dbl_mb(k_singles+i) = 0.0d0
 
152
              dbl_mb(k_doubles+i) = 0.0d0
 
153
            enddo
 
154
 
 
155
           call tce_mrcc_q3vnt2(dbl_mb(k_doubles),d_t2,d_v2,k_t2_offset,
 
156
     1        k_v2_offset,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,2)
 
157
           call tce_mrcc_q3vnt1(dbl_mb(k_singles),d_t1,d_v2,k_t1_offset,
 
158
     1        k_v2_offset,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,2)
 
159
           call tce_mrcc_q3fnt2(dbl_mb(k_singles),d_f1,d_t2,k_f1_offset,
 
160
     +           k_t2_offset,t_h1b,t_h2b,t_h3b,t_p4b,t_p5b,t_p6b,2)
 
161
 
 
162
            if (restricted) then
 
163
              factor = 2.0d0
 
164
            else
 
165
              factor = 1.0d0
 
166
            endif
 
167
            if ((t_p4b.eq.t_p5b).and.(t_p5b.eq.t_p6b)) then
 
168
              factor = factor / 6.0d0
 
169
            else if ((t_p4b.eq.t_p5b).or.(t_p5b.eq.t_p6b)) then
 
170
              factor = factor / 2.0d0
 
171
            endif
 
172
            if ((t_h1b.eq.t_h2b).and.(t_h2b.eq.t_h3b)) then
 
173
              factor = factor / 6.0d0
 
174
            else if ((t_h1b.eq.t_h2b).or.(t_h2b.eq.t_h3b)) then
 
175
              factor = factor / 2.0d0
 
176
            endif
 
177
c
 
178
c factor = [ 1/36, 1/18, 1/12, 1/6, 1/4, 1/3, 1/2, 1, 2]
 
179
c
 
180
            i = 0
 
181
            do t_p4 = 1, range_p4
 
182
             denom_p4 = dbl_mb(offset_p4+t_p4)
 
183
             do t_p5 = 1, range_p5
 
184
              denom_p5 = dbl_mb(offset_p5+t_p5)
 
185
              do t_p6 = 1, range_p6
 
186
               denom_p6 = dbl_mb(offset_p6+t_p6)
 
187
               do t_h1 = 1, range_h1
 
188
                denom_h1 = dbl_mb(offset_h1+t_h1)
 
189
                do t_h2 = 1, range_h2
 
190
                 denom_h2 = dbl_mb(offset_h2+t_h2)
 
191
                 do t_h3 = 1, range_h3
 
192
                  denom_h3 = dbl_mb(offset_h3+t_h3)
 
193
 
 
194
 
 
195
                  if(.not.lusesamefock_nonit) then
 
196
                  denom = 1.0d0 / ((denom_h1 + denom_h2 + denom_h3)
 
197
     +                           - (denom_p4 + denom_p5 + denom_p6))
 
198
 
 
199
 
 
200
                  else
 
201
                  orbspin(1) = int_mb(k_spin+t_p4b-1)-1
 
202
                  orbspin(2) = int_mb(k_spin+t_p5b-1)-1
 
203
                  orbspin(3) = int_mb(k_spin+t_p6b-1)-1
 
204
                  orbspin(4) = int_mb(k_spin+t_h1b-1)-1
 
205
                  orbspin(5) = int_mb(k_spin+t_h2b-1)-1
 
206
                  orbspin(6) = int_mb(k_spin+t_h3b-1)-1
 
207
 
 
208
                  orbindex(1) = (1-orbspin(1)+int_mb(k_mo_indexm(iref)+
 
209
     1           int_mb(k_offset+t_p4b-1)+t_p4-1))/2
 
210
                  orbindex(2) = (1-orbspin(2)+int_mb(k_mo_indexm(iref)+
 
211
     1           int_mb(k_offset+t_p5b-1)+t_p5-1))/2
 
212
                  orbindex(3) = (1-orbspin(3)+int_mb(k_mo_indexm(iref)+
 
213
     1           int_mb(k_offset+t_p6b-1)+t_p6-1))/2
 
214
                  orbindex(4) = (1-orbspin(4)+int_mb(k_mo_indexm(iref)+
 
215
     1           int_mb(k_offset+t_h1b-1)+t_h1-1))/2
 
216
                  orbindex(5) = (1-orbspin(5)+int_mb(k_mo_indexm(iref)+
 
217
     1           int_mb(k_offset+t_h2b-1)+t_h2-1))/2
 
218
                  orbindex(6) = (1-orbspin(6)+int_mb(k_mo_indexm(iref)+
 
219
     1           int_mb(k_offset+t_h3b-1)+t_h3-1))/2
 
220
 
 
221
                 do j=1,6
 
222
                  blck(j) = orbinblck(orbindex(j),orbspin(j)+1,1)
 
223
                  off(j) = offsetinblck(orbindex(j),orbspin(j)+1,1)
 
224
                 enddo
 
225
 
 
226
                denom = 1.0d0/(
 
227
     +                  -dbl_mb(k_evl_sortedm(1)+int_mb(k_offsetm(1)
 
228
     +                         +blck(1)-1)+off(1))
 
229
     +                  -dbl_mb(k_evl_sortedm(1)+int_mb(k_offsetm(1)
 
230
     +                         +blck(2)-1)+off(2))
 
231
     +                  -dbl_mb(k_evl_sortedm(1)+int_mb(k_offsetm(1)
 
232
     +                         +blck(3)-1)+off(3))
 
233
     +                  +dbl_mb(k_evl_sortedm(1)+int_mb(k_offsetm(1)
 
234
     +                         +blck(4)-1)+off(4))
 
235
     +                  +dbl_mb(k_evl_sortedm(1)+int_mb(k_offsetm(1)
 
236
     +                         +blck(5)-1)+off(5))
 
237
     +                  +dbl_mb(k_evl_sortedm(1)+int_mb(k_offsetm(1)
 
238
     +                         +blck(6)-1)+off(6)))
 
239
                endif ! lusesamefock_nonit
 
240
 
 
241
                if((abs(denom)).gt.1.0d2) then
 
242
c                 if(nodezero) 
 
243
                 write(LuOut,*)
 
244
     +           'Warning:Denominator is very low. 1/D=',abs(denom)
 
245
c                 if(nodezero) 
 
246
                 call util_flush(LuOut)
 
247
                endif
 
248
               
 
249
 
 
250
c                  energy1 = energy1 + factor*denom
 
251
c     1                    * dbl_mb(k_doubles+i)*dbl_mb(k_doubles+i)
 
252
c                  energy2 = energy2 + factor*denom*dbl_mb(k_doubles+i)
 
253
c     1                    * (dbl_mb(k_doubles+i)+dbl_mb(k_singles+i))
 
254
                  energy3 = energy3 + factor*denom*dbl_mb(k_doubles+i)
 
255
     +                    * (dbl_mb(k_doubles+i)+dbl_mb(k_singles+i))
 
256
 
 
257
 
 
258
                  i = i + 1
 
259
                 enddo
 
260
                enddo
 
261
               enddo
 
262
              enddo
 
263
             enddo
 
264
            enddo
 
265
c            if (.not.MA_POP_STACK(l_q3fnt2)) 
 
266
c     1        call errquit('ccsd_t',2,MA_ERR)
 
267
            if (.not.MA_POP_STACK(l_doubles)) 
 
268
     1        call errquit('ccsd_t',3,MA_ERR)
 
269
            if (.not.MA_POP_STACK(l_singles)) 
 
270
     1        call errquit('ccsd_t',4,MA_ERR)
 
271
            ENDIF
 
272
            endif
 
273
            endif
 
274
            endif
 
275
             if(lusesub) then
 
276
              next = NXTASKsub(nprocs,1,mypgid)
 
277
             else
 
278
              next = nxtask(nprocs,1)
 
279
             endif
 
280
            endif
 
281
            count = count + 1
 
282
           enddo
 
283
          enddo
 
284
         enddo
 
285
        enddo
 
286
       enddo
 
287
      enddo
 
288
      if(lusesub) then
 
289
      next = NXTASKsub(-nprocs,1,mypgid)
 
290
      else
 
291
      next = nxtask(-nprocs,1)
 
292
      call ga_sync()
 
293
      endif
 
294
c      energy(1) = energy1
 
295
c      energy(2) = energy2
 
296
      energy(3) = energy3
 
297
      if(lusesub) then
 
298
      call ga_pgroup_sync(mypgid)
 
299
c      call ga_pgroup_dgop(p_handle, type, buf, lenbuf, op)
 
300
      call ga_pgroup_dgop(mypgid, mt_dbl,energy,3,'+')
 
301
      else 
 
302
      call ga_dgop(mt_dbl,energy,3,'+')
 
303
      endif
 
304
c      energy1 = energy(1)
 
305
c      energy2 = energy(2)
 
306
      energy3 = energy(3)
 
307
 
 
308
c       if(lusesub)then
 
309
c       if(ga_pgroup_nodeid(mypgid).eq.0) then
 
310
cc      write(*,'(A,F17.10,A,F17.10,A,F17.10)')"energy1 ",energy1,
 
311
cc     +   " energy2 ",energy2," energy3 ",energy3
 
312
c      write(*,'(A,F17.10,A,F17.10,A,F17.10)')"energy3 ",energy3
 
313
c       endif
 
314
c       else
 
315
c      if(nodezero) then
 
316
cc      write(*,'(A,F17.10,A,F17.10,A,F17.10)')"energy1 ",energy1,
 
317
cc     +   " energy2 ",energy2," energy3 ",energy3
 
318
c      write(*,'(A,F17.10,A,F17.10,A,F17.10)')"energy3 ",energy3
 
319
c      endif
 
320
c
 
321
c      endif
 
322
c
 
323
      return
 
324
      end
 
325