1 SUBROUTINE ccsd2_q(d_t1,k_t1_offset,d_t2,k_t2_offset, 2 1 d_v2,k_v2_offset,d_y2,k_y2_offset, 3 2 d_e,k_e_offset,energy) 4C 5C $Id$ 6C 7 IMPLICIT NONE 8#include "global.fh" 9#include "mafdecls.fh" 10#include "util.fh" 11#include "errquit.fh" 12#include "tce.fh" 13#include "tce_main.fh" 14 integer d_t1 15 integer k_t1_offset 16 integer d_t2 17 integer k_t2_offset 18 integer d_v2 19 integer k_v2_offset 20 integer d_y2 21 integer k_y2_offset 22 integer d_e 23 integer k_e_offset 24 integer t_h1b, t_h1 25 integer t_h2b, t_h2 26 integer t_h3b, t_h3 27 integer t_h4b, t_h4 28 integer t_p5b, t_p5 29 integer t_p6b, t_p6 30 integer t_p7b, t_p7 31 integer t_p8b, t_p8 32 integer k_right,l_right 33 integer k_left,l_left 34 integer k_left_sorted,l_left_sorted 35 integer size,i 36 integer g_energy 37 integer nxtask 38 integer next 39 integer nprocs 40 integer count 41 integer d_i1_1,d_i1_2 42 integer k_i1_offset_1,l_i1_offset_1 43 integer k_i1_offset_2,l_i1_offset_2 44 double precision energy 45 double precision factor 46 external nxtask 47c 48c Caution! k_right & k_left are not even allocated yet 49c but they will not be used. 50c 51 call ccsd2_q_right(dbl_mb(k_right),d_i1_1,d_i1_2, 52 1 d_t1,d_t2,d_v2,k_i1_offset_1,k_i1_offset_2, 53 2 k_t1_offset,k_t2_offset,k_v2_offset,l_i1_offset_1, 54 3 l_i1_offset_2,t_h1b,t_h2b,t_h3b,t_h4b, 55 4 t_p5b,t_p6b,t_p7b,t_p8b,1) 56 call ccsd2_q_left(dbl_mb(k_left),d_v2,d_y2,k_v2_offset, 57 1 k_y2_offset,t_h1b,t_h2b,t_h3b,t_h4b, 58 2 t_p5b,t_p6b,t_p7b,t_p8b,1) 59c 60c Accumulate 61c 62 if (.not.ga_create(mt_dbl,1,1,'perturbative',1,1,g_energy)) 63 1 call errquit('ccsd_q: GA problem',0,GA_ERR) 64 nprocs = GA_NNODES() 65 count = 0 66 next = nxtask(nprocs,1) 67 energy=0.0d0 68 do t_p5b = noab+1,noab+nvab 69 do t_p6b = t_p5b,noab+nvab 70 do t_p7b = t_p6b,noab+nvab 71 do t_p8b = t_p7b,noab+nvab 72 do t_h1b = 1,noab 73 do t_h2b = t_h1b,noab 74 do t_h3b = t_h2b,noab 75 do t_h4b = t_h3b,noab 76 if (next.eq.count) then 77 if (int_mb(k_spin+t_p5b-1) 78 1 +int_mb(k_spin+t_p6b-1) 79 2 +int_mb(k_spin+t_p7b-1) 80 3 +int_mb(k_spin+t_p8b-1) 81 4 .eq.int_mb(k_spin+t_h1b-1) 82 5 +int_mb(k_spin+t_h2b-1) 83 6 +int_mb(k_spin+t_h3b-1) 84 7 +int_mb(k_spin+t_h4b-1)) then 85 if ((.not.restricted).or. 86 1 (int_mb(k_spin+t_p5b-1) 87 2 +int_mb(k_spin+t_p6b-1) 88 3 +int_mb(k_spin+t_p7b-1) 89 4 +int_mb(k_spin+t_p8b-1) 90 5 +int_mb(k_spin+t_h1b-1) 91 6 +int_mb(k_spin+t_h2b-1) 92 7 +int_mb(k_spin+t_h3b-1) 93 8 +int_mb(k_spin+t_h4b-1).le.12)) then 94 if (ieor(int_mb(k_sym+t_p5b-1), 95 1 ieor(int_mb(k_sym+t_p6b-1), 96 2 ieor(int_mb(k_sym+t_p7b-1), 97 3 ieor(int_mb(k_sym+t_p8b-1), 98 4 ieor(int_mb(k_sym+t_h1b-1), 99 5 ieor(int_mb(k_sym+t_h2b-1), 100 6 ieor(int_mb(k_sym+t_h3b-1), 101 7 int_mb(k_sym+t_h4b-1)))))))).eq.0) then 102 size = int_mb(k_range+t_p5b-1) 103 1 * int_mb(k_range+t_p6b-1) 104 2 * int_mb(k_range+t_p7b-1) 105 3 * int_mb(k_range+t_p8b-1) 106 4 * int_mb(k_range+t_h1b-1) 107 5 * int_mb(k_range+t_h2b-1) 108 6 * int_mb(k_range+t_h3b-1) 109 7 * int_mb(k_range+t_h4b-1) 110 if (.not.MA_PUSH_GET(mt_dbl,size,'right moment 2,3', 111 1 l_right,k_right)) call errquit('ccsd_q',3,MA_ERR) 112 if (.not.MA_PUSH_GET(mt_dbl,size,'left moment 2,3', 113 1 l_left,k_left)) call errquit('ccsd_q',3,MA_ERR) 114 if (.not.MA_PUSH_GET(mt_dbl,size,'left moment 2,3', 115 1 l_left_sorted,k_left_sorted)) call errquit('ccsd_q',3, 116 2 MA_ERR) 117 do i = 1, size 118 dbl_mb(k_right+i-1) = 0.0d0 119 dbl_mb(k_left+i-1) = 0.0d0 120 enddo 121 call ccsd2_q_right(dbl_mb(k_right),d_i1_1,d_i1_2, 122 1 d_t1,d_t2,d_v2,k_i1_offset_1,k_i1_offset_2, 123 2 k_t1_offset,k_t2_offset,k_v2_offset,l_i1_offset_1, 124 3 l_i1_offset_2,t_h1b,t_h2b,t_h3b,t_h4b, 125 4 t_p5b,t_p6b,t_p7b,t_p8b,2) 126 call ccsd2_q_left(dbl_mb(k_left),d_v2,d_y2,k_v2_offset, 127 1 k_y2_offset,t_h1b,t_h2b,t_h3b,t_h4b, 128 2 t_p5b,t_p6b,t_p7b,t_p8b,2) 129 call tce_sort_8(dbl_mb(k_left),dbl_mb(k_left_sorted), 130 1 int_mb(k_range+t_h1b-1),int_mb(k_range+t_h2b-1), 131 2 int_mb(k_range+t_h3b-1),int_mb(k_range+t_h4b-1), 132 3 int_mb(k_range+t_p5b-1),int_mb(k_range+t_p6b-1), 133 4 int_mb(k_range+t_p7b-1),int_mb(k_range+t_p8b-1), 134 5 5,6,7,8,1,2,3,4,1.0d0) 135 if ((restricted).and. 136 1 (int_mb(k_spin+t_p5b-1) 137 2 +int_mb(k_spin+t_p6b-1) 138 3 +int_mb(k_spin+t_p7b-1) 139 4 +int_mb(k_spin+t_p8b-1) 140 5 +int_mb(k_spin+t_h1b-1) 141 6 +int_mb(k_spin+t_h2b-1) 142 7 +int_mb(k_spin+t_h3b-1) 143 8 +int_mb(k_spin+t_h4b-1).ne.12)) then 144 factor = 2.0d0 145 else 146 factor = 1.0d0 147 endif 148 if ((t_p5b.eq.t_p6b).and.(t_p6b.eq.t_p7b) 149 1 .and.(t_p7b.eq.t_p8b)) then 150 factor = factor / 24.0d0 151 else if ((t_p5b.eq.t_p6b).and.(t_p6b.eq.t_p7b)) then 152 factor = factor / 6.0d0 153 else if ((t_p6b.eq.t_p7b).and.(t_p7b.eq.t_p8b)) then 154 factor = factor / 6.0d0 155 else if ((t_p5b.eq.t_p6b).and.(t_p7b.eq.t_p8b)) then 156 factor = factor / 4.0d0 157 else if (t_p5b.eq.t_p6b) then 158 factor = factor / 2.0d0 159 else if (t_p6b.eq.t_p7b) then 160 factor = factor / 2.0d0 161 else if (t_p7b.eq.t_p8b) then 162 factor = factor / 2.0d0 163 endif 164 if ((t_h1b.eq.t_h2b).and.(t_h2b.eq.t_h3b) 165 1 .and.(t_h3b.eq.t_h4b)) then 166 factor = factor / 24.0d0 167 else if ((t_h1b.eq.t_h2b).and.(t_h2b.eq.t_h3b)) then 168 factor = factor / 6.0d0 169 else if ((t_h2b.eq.t_h3b).and.(t_h3b.eq.t_h4b)) then 170 factor = factor / 6.0d0 171 else if ((t_h1b.eq.t_h2b).and.(t_h3b.eq.t_h4b)) then 172 factor = factor / 4.0d0 173 else if (t_h1b.eq.t_h2b) then 174 factor = factor / 2.0d0 175 else if (t_h2b.eq.t_h3b) then 176 factor = factor / 2.0d0 177 else if (t_h3b.eq.t_h4b) then 178 factor = factor / 2.0d0 179 endif 180 i = 0 181 do t_p5 = 1, int_mb(k_range+t_p5b-1) 182 do t_p6 = 1, int_mb(k_range+t_p6b-1) 183 do t_p7 = 1, int_mb(k_range+t_p7b-1) 184 do t_p8 = 1, int_mb(k_range+t_p8b-1) 185 do t_h1 = 1, int_mb(k_range+t_h1b-1) 186 do t_h2 = 1, int_mb(k_range+t_h2b-1) 187 do t_h3 = 1, int_mb(k_range+t_h3b-1) 188 do t_h4 = 1, int_mb(k_range+t_h4b-1) 189 i = i + 1 190 energy = energy + factor * dbl_mb(k_right+i-1) 191 1 * dbl_mb(k_left_sorted+i-1) 192 2 / (-dbl_mb(k_evl_sorted+int_mb(k_offset+t_p5b-1)+t_p5-1) 193 3 -dbl_mb(k_evl_sorted+int_mb(k_offset+t_p6b-1)+t_p6-1) 194 4 -dbl_mb(k_evl_sorted+int_mb(k_offset+t_p7b-1)+t_p7-1) 195 5 -dbl_mb(k_evl_sorted+int_mb(k_offset+t_p8b-1)+t_p8-1) 196 6 +dbl_mb(k_evl_sorted+int_mb(k_offset+t_h1b-1)+t_h1-1) 197 7 +dbl_mb(k_evl_sorted+int_mb(k_offset+t_h2b-1)+t_h2-1) 198 8 +dbl_mb(k_evl_sorted+int_mb(k_offset+t_h3b-1)+t_h3-1) 199 9 +dbl_mb(k_evl_sorted+int_mb(k_offset+t_h4b-1)+t_h4-1)) 200 enddo 201 enddo 202 enddo 203 enddo 204 enddo 205 enddo 206 enddo 207 enddo 208 if (.not.MA_POP_STACK(l_left_sorted)) 209 1 call errquit('ccsd_q',6,MA_ERR) 210 if (.not.MA_POP_STACK(l_left)) 211 1 call errquit('ccsd_q',6,MA_ERR) 212 if (.not.MA_POP_STACK(l_right)) 213 1 call errquit('ccsd_q',6,MA_ERR) 214 endif 215 endif 216 endif 217 next = nxtask(nprocs,1) 218 endif 219 count = count + 1 220 enddo 221 enddo 222 enddo 223 enddo 224 enddo 225 enddo 226 enddo 227 enddo 228 next = nxtask(-nprocs,1) 229 call ccsd2_q_left(dbl_mb(k_left),d_v2,d_y2,k_v2_offset, 230 1 k_y2_offset,t_h1b,t_h2b,t_h3b,t_h4b, 231 2 t_p5b,t_p6b,t_p7b,t_p8b,3) 232 call ccsd2_q_right(dbl_mb(k_right),d_i1_1,d_i1_2, 233 1 d_t1,d_t2,d_v2,k_i1_offset_1,k_i1_offset_2, 234 2 k_t1_offset,k_t2_offset,k_v2_offset,l_i1_offset_1, 235 3 l_i1_offset_2,t_h1b,t_h2b,t_h3b,t_h4b, 236 4 t_p5b,t_p6b,t_p7b,t_p8b,3) 237 call ga_zero(g_energy) 238 call ga_acc(g_energy,1,1,1,1,energy,1,1.0d0) 239 call ga_sync() 240 call ga_get(g_energy,1,1,1,1,energy,1) 241 if (.not.ga_destroy(g_energy)) 242 1 call errquit('ccsd_q: GA problem',1,GA_ERR) 243 return 244 end 245