1 SUBROUTINE ccsdt_o1(d_i0,d_o1,d_t1,d_t2,k_i0_offset,k_o1_offset,k_ 2 &t1_offset,k_t2_offset) 3C $Id$ 4C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 5C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 6C i0 ( p2 h1 )_o + = 1 * o ( p2 h1 )_o 7C i0 ( p2 h1 )_to + = -1 * Sum ( h5 ) * t ( p2 h5 )_t * i1 ( h5 h1 )_o 8C i1 ( h5 h1 )_o + = 1 * o ( h5 h1 )_o 9C i1 ( h5 h1 )_ot + = 1 * Sum ( p3 ) * o ( h5 p3 )_o * t ( p3 h1 )_t 10C i0 ( p2 h1 )_to + = 1 * Sum ( p3 ) * o ( p2 p3 )_o * t ( p3 h1 )_t 11C i0 ( p2 h1 )_to + = 1 * Sum ( p4 h3 ) * o ( h3 p4 )_o * t ( p2 p4 h1 h3 )_t 12 IMPLICIT NONE 13#include "global.fh" 14#include "mafdecls.fh" 15#include "util.fh" 16#include "errquit.fh" 17#include "tce.fh" 18 INTEGER d_i0 19 INTEGER k_i0_offset 20 INTEGER d_o1 21 INTEGER k_o1_offset 22 INTEGER d_t1 23 INTEGER k_t1_offset 24 INTEGER d_i1 25 INTEGER k_i1_offset 26 INTEGER d_t2 27 INTEGER k_t2_offset 28 INTEGER l_i1_offset 29 INTEGER size_i1 30 CHARACTER*255 filename 31 CALL ccsdt_o1_1(d_o1,k_o1_offset,d_i0,k_i0_offset) 32 CALL OFFSET_ccsdt_o1_2_1(l_i1_offset,k_i1_offset,size_i1) 33 CALL TCE_FILENAME('ccsdt_o1_2_1_i1',filename) 34 CALL CREATEFILE(filename,d_i1,size_i1) 35 CALL ccsdt_o1_2_1(d_o1,k_o1_offset,d_i1,k_i1_offset) 36 CALL ccsdt_o1_2_2(d_o1,k_o1_offset,d_t1,k_t1_offset,d_i1,k_i1_offs 37 &et) 38 CALL RECONCILEFILE(d_i1,size_i1) 39 CALL ccsdt_o1_2(d_t1,k_t1_offset,d_i1,k_i1_offset,d_i0,k_i0_offset 40 &) 41 CALL DELETEFILE(d_i1) 42 IF (.not.MA_POP_STACK(l_i1_offset)) CALL ERRQUIT('ccsdt_o1',-1,MA_ 43 &ERR) 44 CALL ccsdt_o1_3(d_o1,k_o1_offset,d_t1,k_t1_offset,d_i0,k_i0_offset 45 &) 46 CALL ccsdt_o1_4(d_o1,k_o1_offset,d_t2,k_t2_offset,d_i0,k_i0_offset 47 &) 48 RETURN 49 END 50 SUBROUTINE ccsdt_o1_1(d_a,k_a_offset,d_c,k_c_offset) 51C $Id$ 52C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 53C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 54C i0 ( p2 h1 )_o + = 1 * o ( p2 h1 )_o 55 IMPLICIT NONE 56#include "global.fh" 57#include "mafdecls.fh" 58#include "sym.fh" 59#include "errquit.fh" 60#include "tce.fh" 61 INTEGER d_a 62 INTEGER k_a_offset 63 INTEGER d_c 64 INTEGER k_c_offset 65 INTEGER NXTASK 66 INTEGER next 67 INTEGER nprocs 68 INTEGER count 69 INTEGER p2b 70 INTEGER h1b 71 INTEGER dimc 72 INTEGER p2b_1 73 INTEGER h1b_1 74 INTEGER dim_common 75 INTEGER dima_sort 76 INTEGER dima 77 INTEGER l_a_sort 78 INTEGER k_a_sort 79 INTEGER l_a 80 INTEGER k_a 81 INTEGER l_c 82 INTEGER k_c 83 EXTERNAL NXTASK 84 nprocs = GA_NNODES() 85 count = 0 86 next = NXTASK(nprocs,1) 87 DO p2b = noab+1,noab+nvab 88 DO h1b = 1,noab 89 IF (next.eq.count) THEN 90 IF ((.not.restricted).or.(int_mb(k_spin+p2b-1)+int_mb(k_spin+h1b-1 91 &).ne.4)) THEN 92 IF (int_mb(k_spin+p2b-1) .eq. int_mb(k_spin+h1b-1)) THEN 93 IF (ieor(int_mb(k_sym+p2b-1),int_mb(k_sym+h1b-1)) .eq. irrep_o) TH 94 &EN 95 dimc = int_mb(k_range+p2b-1) * int_mb(k_range+h1b-1) 96 CALL TCE_RESTRICTED_2(p2b,h1b,p2b_1,h1b_1) 97 dim_common = 1 98 dima_sort = int_mb(k_range+p2b-1) * int_mb(k_range+h1b-1) 99 dima = dim_common * dima_sort 100 IF (dima .gt. 0) THEN 101 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a_sort,k_a_sort)) CALL 102 & ERRQUIT('ccsdt_o1_1',0,MA_ERR) 103 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a,k_a)) CALL ERRQUIT(' 104 &ccsdt_o1_1',1,MA_ERR) 105 CALL GET_HASH_BLOCK(d_a,dbl_mb(k_a),dima,int_mb(k_a_offset),(h1b_1 106 & - 1 + (noab+nvab) * (p2b_1 - 1))) 107 CALL TCE_SORT_2(dbl_mb(k_a),dbl_mb(k_a_sort),int_mb(k_range+p2b-1) 108 &,int_mb(k_range+h1b-1),2,1,1.0d0) 109 IF (.not.MA_POP_STACK(l_a)) CALL ERRQUIT('ccsdt_o1_1',2,MA_ERR) 110 IF (.not.MA_PUSH_GET(mt_dbl,dimc,'noname',l_c,k_c)) CALL ERRQUIT(' 111 &ccsdt_o1_1',3,MA_ERR) 112 CALL TCE_SORT_2(dbl_mb(k_a_sort),dbl_mb(k_c),int_mb(k_range+h1b-1) 113 &,int_mb(k_range+p2b-1),2,1,1.0d0) 114 CALL ADD_HASH_BLOCK(d_c,dbl_mb(k_c),dimc,int_mb(k_c_offset),(h1b - 115 & 1 + noab * (p2b - noab - 1))) 116 IF (.not.MA_POP_STACK(l_c)) CALL ERRQUIT('ccsdt_o1_1',4,MA_ERR) 117 IF (.not.MA_POP_STACK(l_a_sort)) CALL ERRQUIT('ccsdt_o1_1',5,MA_ER 118 &R) 119 END IF 120 END IF 121 END IF 122 END IF 123 next = NXTASK(nprocs,1) 124 END IF 125 count = count + 1 126 END DO 127 END DO 128 next = NXTASK(-nprocs,1) 129 call GA_SYNC() 130 RETURN 131 END 132 SUBROUTINE ccsdt_o1_2(d_a,k_a_offset,d_b,k_b_offset,d_c,k_c_offset 133 &) 134C $Id$ 135C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 136C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 137C i0 ( p2 h1 )_to + = -1 * Sum ( h5 ) * t ( p2 h5 )_t * i1 ( h5 h1 )_o 138 IMPLICIT NONE 139#include "global.fh" 140#include "mafdecls.fh" 141#include "sym.fh" 142#include "errquit.fh" 143#include "tce.fh" 144 INTEGER d_a 145 INTEGER k_a_offset 146 INTEGER d_b 147 INTEGER k_b_offset 148 INTEGER d_c 149 INTEGER k_c_offset 150 INTEGER NXTASK 151 INTEGER next 152 INTEGER nprocs 153 INTEGER count 154 INTEGER p2b 155 INTEGER h1b 156 INTEGER dimc 157 INTEGER l_c_sort 158 INTEGER k_c_sort 159 INTEGER h5b 160 INTEGER p2b_1 161 INTEGER h5b_1 162 INTEGER h5b_2 163 INTEGER h1b_2 164 INTEGER dim_common 165 INTEGER dima_sort 166 INTEGER dima 167 INTEGER dimb_sort 168 INTEGER dimb 169 INTEGER l_a_sort 170 INTEGER k_a_sort 171 INTEGER l_a 172 INTEGER k_a 173 INTEGER l_b_sort 174 INTEGER k_b_sort 175 INTEGER l_b 176 INTEGER k_b 177 INTEGER l_c 178 INTEGER k_c 179 EXTERNAL NXTASK 180 nprocs = GA_NNODES() 181 count = 0 182 next = NXTASK(nprocs,1) 183 DO p2b = noab+1,noab+nvab 184 DO h1b = 1,noab 185 IF (next.eq.count) THEN 186 IF ((.not.restricted).or.(int_mb(k_spin+p2b-1)+int_mb(k_spin+h1b-1 187 &).ne.4)) THEN 188 IF (int_mb(k_spin+p2b-1) .eq. int_mb(k_spin+h1b-1)) THEN 189 IF (ieor(int_mb(k_sym+p2b-1),int_mb(k_sym+h1b-1)) .eq. ieor(irrep_ 190 &t,irrep_o)) THEN 191 dimc = int_mb(k_range+p2b-1) * int_mb(k_range+h1b-1) 192 IF (.not.MA_PUSH_GET(mt_dbl,dimc,'noname',l_c_sort,k_c_sort)) CALL 193 & ERRQUIT('ccsdt_o1_2',0,MA_ERR) 194 CALL DFILL(dimc,0.0d0,dbl_mb(k_c_sort),1) 195 DO h5b = 1,noab 196 IF (int_mb(k_spin+p2b-1) .eq. int_mb(k_spin+h5b-1)) THEN 197 IF (ieor(int_mb(k_sym+p2b-1),int_mb(k_sym+h5b-1)) .eq. irrep_t) TH 198 &EN 199 CALL TCE_RESTRICTED_2(p2b,h5b,p2b_1,h5b_1) 200 CALL TCE_RESTRICTED_2(h5b,h1b,h5b_2,h1b_2) 201 dim_common = int_mb(k_range+h5b-1) 202 dima_sort = int_mb(k_range+p2b-1) 203 dima = dim_common * dima_sort 204 dimb_sort = int_mb(k_range+h1b-1) 205 dimb = dim_common * dimb_sort 206 IF ((dima .gt. 0) .and. (dimb .gt. 0)) THEN 207 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a_sort,k_a_sort)) CALL 208 & ERRQUIT('ccsdt_o1_2',1,MA_ERR) 209 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a,k_a)) CALL ERRQUIT(' 210 &ccsdt_o1_2',2,MA_ERR) 211 CALL GET_HASH_BLOCK(d_a,dbl_mb(k_a),dima,int_mb(k_a_offset),(h5b_1 212 & - 1 + noab * (p2b_1 - noab - 1))) 213 CALL TCE_SORT_2(dbl_mb(k_a),dbl_mb(k_a_sort),int_mb(k_range+p2b-1) 214 &,int_mb(k_range+h5b-1),1,2,1.0d0) 215 IF (.not.MA_POP_STACK(l_a)) CALL ERRQUIT('ccsdt_o1_2',3,MA_ERR) 216 IF (.not.MA_PUSH_GET(mt_dbl,dimb,'noname',l_b_sort,k_b_sort)) CALL 217 & ERRQUIT('ccsdt_o1_2',4,MA_ERR) 218 IF (.not.MA_PUSH_GET(mt_dbl,dimb,'noname',l_b,k_b)) CALL ERRQUIT(' 219 &ccsdt_o1_2',5,MA_ERR) 220 CALL GET_HASH_BLOCK(d_b,dbl_mb(k_b),dimb,int_mb(k_b_offset),(h1b_2 221 & - 1 + noab * (h5b_2 - 1))) 222 CALL TCE_SORT_2(dbl_mb(k_b),dbl_mb(k_b_sort),int_mb(k_range+h5b-1) 223 &,int_mb(k_range+h1b-1),2,1,1.0d0) 224 IF (.not.MA_POP_STACK(l_b)) CALL ERRQUIT('ccsdt_o1_2',6,MA_ERR) 225 CALL DGEMM('T','N',dima_sort,dimb_sort,dim_common,1.0d0,dbl_mb(k_a 226 &_sort),dim_common,dbl_mb(k_b_sort),dim_common,1.0d0,dbl_mb(k_c_sor 227 &t),dima_sort) 228 IF (.not.MA_POP_STACK(l_b_sort)) CALL ERRQUIT('ccsdt_o1_2',7,MA_ER 229 &R) 230 IF (.not.MA_POP_STACK(l_a_sort)) CALL ERRQUIT('ccsdt_o1_2',8,MA_ER 231 &R) 232 END IF 233 END IF 234 END IF 235 END DO 236 IF (.not.MA_PUSH_GET(mt_dbl,dimc,'noname',l_c,k_c)) CALL ERRQUIT(' 237 &ccsdt_o1_2',9,MA_ERR) 238 CALL TCE_SORT_2(dbl_mb(k_c_sort),dbl_mb(k_c),int_mb(k_range+h1b-1) 239 &,int_mb(k_range+p2b-1),2,1,-1.0d0) 240 CALL ADD_HASH_BLOCK(d_c,dbl_mb(k_c),dimc,int_mb(k_c_offset),(h1b - 241 & 1 + noab * (p2b - noab - 1))) 242 IF (.not.MA_POP_STACK(l_c)) CALL ERRQUIT('ccsdt_o1_2',10,MA_ERR) 243 IF (.not.MA_POP_STACK(l_c_sort)) CALL ERRQUIT('ccsdt_o1_2',11,MA_E 244 &RR) 245 END IF 246 END IF 247 END IF 248 next = NXTASK(nprocs,1) 249 END IF 250 count = count + 1 251 END DO 252 END DO 253 next = NXTASK(-nprocs,1) 254 call GA_SYNC() 255 RETURN 256 END 257 SUBROUTINE ccsdt_o1_2_1(d_a,k_a_offset,d_c,k_c_offset) 258C $Id$ 259C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 260C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 261C i1 ( h5 h1 )_o + = 1 * o ( h5 h1 )_o 262 IMPLICIT NONE 263#include "global.fh" 264#include "mafdecls.fh" 265#include "sym.fh" 266#include "errquit.fh" 267#include "tce.fh" 268 INTEGER d_a 269 INTEGER k_a_offset 270 INTEGER d_c 271 INTEGER k_c_offset 272 INTEGER NXTASK 273 INTEGER next 274 INTEGER nprocs 275 INTEGER count 276 INTEGER h5b 277 INTEGER h1b 278 INTEGER dimc 279 INTEGER h5b_1 280 INTEGER h1b_1 281 INTEGER dim_common 282 INTEGER dima_sort 283 INTEGER dima 284 INTEGER l_a_sort 285 INTEGER k_a_sort 286 INTEGER l_a 287 INTEGER k_a 288 INTEGER l_c 289 INTEGER k_c 290 EXTERNAL NXTASK 291 nprocs = GA_NNODES() 292 count = 0 293 next = NXTASK(nprocs,1) 294 DO h5b = 1,noab 295 DO h1b = 1,noab 296 IF (next.eq.count) THEN 297 IF ((.not.restricted).or.(int_mb(k_spin+h5b-1)+int_mb(k_spin+h1b-1 298 &).ne.4)) THEN 299 IF (int_mb(k_spin+h5b-1) .eq. int_mb(k_spin+h1b-1)) THEN 300 IF (ieor(int_mb(k_sym+h5b-1),int_mb(k_sym+h1b-1)) .eq. irrep_o) TH 301 &EN 302 dimc = int_mb(k_range+h5b-1) * int_mb(k_range+h1b-1) 303 CALL TCE_RESTRICTED_2(h5b,h1b,h5b_1,h1b_1) 304 dim_common = 1 305 dima_sort = int_mb(k_range+h5b-1) * int_mb(k_range+h1b-1) 306 dima = dim_common * dima_sort 307 IF (dima .gt. 0) THEN 308 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a_sort,k_a_sort)) CALL 309 & ERRQUIT('ccsdt_o1_2_1',0,MA_ERR) 310 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a,k_a)) CALL ERRQUIT(' 311 &ccsdt_o1_2_1',1,MA_ERR) 312 CALL GET_HASH_BLOCK(d_a,dbl_mb(k_a),dima,int_mb(k_a_offset),(h1b_1 313 & - 1 + (noab+nvab) * (h5b_1 - 1))) 314 CALL TCE_SORT_2(dbl_mb(k_a),dbl_mb(k_a_sort),int_mb(k_range+h5b-1) 315 &,int_mb(k_range+h1b-1),2,1,1.0d0) 316 IF (.not.MA_POP_STACK(l_a)) CALL ERRQUIT('ccsdt_o1_2_1',2,MA_ERR) 317 IF (.not.MA_PUSH_GET(mt_dbl,dimc,'noname',l_c,k_c)) CALL ERRQUIT(' 318 &ccsdt_o1_2_1',3,MA_ERR) 319 CALL TCE_SORT_2(dbl_mb(k_a_sort),dbl_mb(k_c),int_mb(k_range+h1b-1) 320 &,int_mb(k_range+h5b-1),2,1,1.0d0) 321 CALL ADD_HASH_BLOCK(d_c,dbl_mb(k_c),dimc,int_mb(k_c_offset),(h1b - 322 & 1 + noab * (h5b - 1))) 323 IF (.not.MA_POP_STACK(l_c)) CALL ERRQUIT('ccsdt_o1_2_1',4,MA_ERR) 324 IF (.not.MA_POP_STACK(l_a_sort)) CALL ERRQUIT('ccsdt_o1_2_1',5,MA_ 325 &ERR) 326 END IF 327 END IF 328 END IF 329 END IF 330 next = NXTASK(nprocs,1) 331 END IF 332 count = count + 1 333 END DO 334 END DO 335 next = NXTASK(-nprocs,1) 336 call GA_SYNC() 337 RETURN 338 END 339 SUBROUTINE OFFSET_ccsdt_o1_2_1(l_a_offset,k_a_offset,size) 340C $Id$ 341C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 342C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 343C i1 ( h5 h1 )_o 344 IMPLICIT NONE 345#include "global.fh" 346#include "mafdecls.fh" 347#include "sym.fh" 348#include "errquit.fh" 349#include "tce.fh" 350 INTEGER l_a_offset 351 INTEGER k_a_offset 352 INTEGER size 353 INTEGER length 354 INTEGER addr 355 INTEGER h5b 356 INTEGER h1b 357 length = 0 358 DO h5b = 1,noab 359 DO h1b = 1,noab 360 IF (int_mb(k_spin+h5b-1) .eq. int_mb(k_spin+h1b-1)) THEN 361 IF (ieor(int_mb(k_sym+h5b-1),int_mb(k_sym+h1b-1)) .eq. irrep_o) TH 362 &EN 363 IF ((.not.restricted).or.(int_mb(k_spin+h5b-1)+int_mb(k_spin+h1b-1 364 &).ne.4)) THEN 365 length = length + 1 366 END IF 367 END IF 368 END IF 369 END DO 370 END DO 371 IF (.not.MA_PUSH_GET(mt_int,2*length+1,'noname',l_a_offset,k_a_off 372 &set)) CALL ERRQUIT('ccsdt_o1_2_1',0,MA_ERR) 373 int_mb(k_a_offset) = length 374 addr = 0 375 size = 0 376 DO h5b = 1,noab 377 DO h1b = 1,noab 378 IF (int_mb(k_spin+h5b-1) .eq. int_mb(k_spin+h1b-1)) THEN 379 IF (ieor(int_mb(k_sym+h5b-1),int_mb(k_sym+h1b-1)) .eq. irrep_o) TH 380 &EN 381 IF ((.not.restricted).or.(int_mb(k_spin+h5b-1)+int_mb(k_spin+h1b-1 382 &).ne.4)) THEN 383 addr = addr + 1 384 int_mb(k_a_offset+addr) = h1b - 1 + noab * (h5b - 1) 385 int_mb(k_a_offset+length+addr) = size 386 size = size + int_mb(k_range+h5b-1) * int_mb(k_range+h1b-1) 387 END IF 388 END IF 389 END IF 390 END DO 391 END DO 392 RETURN 393 END 394 SUBROUTINE ccsdt_o1_2_2(d_a,k_a_offset,d_b,k_b_offset,d_c,k_c_offs 395 &et) 396C $Id$ 397C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 398C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 399C i1 ( h5 h1 )_ot + = 1 * Sum ( p3 ) * o ( h5 p3 )_o * t ( p3 h1 )_t 400 IMPLICIT NONE 401#include "global.fh" 402#include "mafdecls.fh" 403#include "sym.fh" 404#include "errquit.fh" 405#include "tce.fh" 406 INTEGER d_a 407 INTEGER k_a_offset 408 INTEGER d_b 409 INTEGER k_b_offset 410 INTEGER d_c 411 INTEGER k_c_offset 412 INTEGER NXTASK 413 INTEGER next 414 INTEGER nprocs 415 INTEGER count 416 INTEGER h5b 417 INTEGER h1b 418 INTEGER dimc 419 INTEGER l_c_sort 420 INTEGER k_c_sort 421 INTEGER p3b 422 INTEGER h5b_1 423 INTEGER p3b_1 424 INTEGER p3b_2 425 INTEGER h1b_2 426 INTEGER dim_common 427 INTEGER dima_sort 428 INTEGER dima 429 INTEGER dimb_sort 430 INTEGER dimb 431 INTEGER l_a_sort 432 INTEGER k_a_sort 433 INTEGER l_a 434 INTEGER k_a 435 INTEGER l_b_sort 436 INTEGER k_b_sort 437 INTEGER l_b 438 INTEGER k_b 439 INTEGER l_c 440 INTEGER k_c 441 EXTERNAL NXTASK 442 nprocs = GA_NNODES() 443 count = 0 444 next = NXTASK(nprocs,1) 445 DO h5b = 1,noab 446 DO h1b = 1,noab 447 IF (next.eq.count) THEN 448 IF ((.not.restricted).or.(int_mb(k_spin+h5b-1)+int_mb(k_spin+h1b-1 449 &).ne.4)) THEN 450 IF (int_mb(k_spin+h5b-1) .eq. int_mb(k_spin+h1b-1)) THEN 451 IF (ieor(int_mb(k_sym+h5b-1),int_mb(k_sym+h1b-1)) .eq. ieor(irrep_ 452 &o,irrep_t)) THEN 453 dimc = int_mb(k_range+h5b-1) * int_mb(k_range+h1b-1) 454 IF (.not.MA_PUSH_GET(mt_dbl,dimc,'noname',l_c_sort,k_c_sort)) CALL 455 & ERRQUIT('ccsdt_o1_2_2',0,MA_ERR) 456 CALL DFILL(dimc,0.0d0,dbl_mb(k_c_sort),1) 457 DO p3b = noab+1,noab+nvab 458 IF (int_mb(k_spin+h5b-1) .eq. int_mb(k_spin+p3b-1)) THEN 459 IF (ieor(int_mb(k_sym+h5b-1),int_mb(k_sym+p3b-1)) .eq. irrep_o) TH 460 &EN 461 CALL TCE_RESTRICTED_2(h5b,p3b,h5b_1,p3b_1) 462 CALL TCE_RESTRICTED_2(p3b,h1b,p3b_2,h1b_2) 463 dim_common = int_mb(k_range+p3b-1) 464 dima_sort = int_mb(k_range+h5b-1) 465 dima = dim_common * dima_sort 466 dimb_sort = int_mb(k_range+h1b-1) 467 dimb = dim_common * dimb_sort 468 IF ((dima .gt. 0) .and. (dimb .gt. 0)) THEN 469 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a_sort,k_a_sort)) CALL 470 & ERRQUIT('ccsdt_o1_2_2',1,MA_ERR) 471 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a,k_a)) CALL ERRQUIT(' 472 &ccsdt_o1_2_2',2,MA_ERR) 473 CALL GET_HASH_BLOCK(d_a,dbl_mb(k_a),dima,int_mb(k_a_offset),(p3b_1 474 & - 1 + (noab+nvab) * (h5b_1 - 1))) 475 CALL TCE_SORT_2(dbl_mb(k_a),dbl_mb(k_a_sort),int_mb(k_range+h5b-1) 476 &,int_mb(k_range+p3b-1),1,2,1.0d0) 477 IF (.not.MA_POP_STACK(l_a)) CALL ERRQUIT('ccsdt_o1_2_2',3,MA_ERR) 478 IF (.not.MA_PUSH_GET(mt_dbl,dimb,'noname',l_b_sort,k_b_sort)) CALL 479 & ERRQUIT('ccsdt_o1_2_2',4,MA_ERR) 480 IF (.not.MA_PUSH_GET(mt_dbl,dimb,'noname',l_b,k_b)) CALL ERRQUIT(' 481 &ccsdt_o1_2_2',5,MA_ERR) 482 CALL GET_HASH_BLOCK(d_b,dbl_mb(k_b),dimb,int_mb(k_b_offset),(h1b_2 483 & - 1 + noab * (p3b_2 - noab - 1))) 484 CALL TCE_SORT_2(dbl_mb(k_b),dbl_mb(k_b_sort),int_mb(k_range+p3b-1) 485 &,int_mb(k_range+h1b-1),2,1,1.0d0) 486 IF (.not.MA_POP_STACK(l_b)) CALL ERRQUIT('ccsdt_o1_2_2',6,MA_ERR) 487 CALL DGEMM('T','N',dima_sort,dimb_sort,dim_common,1.0d0,dbl_mb(k_a 488 &_sort),dim_common,dbl_mb(k_b_sort),dim_common,1.0d0,dbl_mb(k_c_sor 489 &t),dima_sort) 490 IF (.not.MA_POP_STACK(l_b_sort)) CALL ERRQUIT('ccsdt_o1_2_2',7,MA_ 491 &ERR) 492 IF (.not.MA_POP_STACK(l_a_sort)) CALL ERRQUIT('ccsdt_o1_2_2',8,MA_ 493 &ERR) 494 END IF 495 END IF 496 END IF 497 END DO 498 IF (.not.MA_PUSH_GET(mt_dbl,dimc,'noname',l_c,k_c)) CALL ERRQUIT(' 499 &ccsdt_o1_2_2',9,MA_ERR) 500 CALL TCE_SORT_2(dbl_mb(k_c_sort),dbl_mb(k_c),int_mb(k_range+h1b-1) 501 &,int_mb(k_range+h5b-1),2,1,1.0d0) 502 CALL ADD_HASH_BLOCK(d_c,dbl_mb(k_c),dimc,int_mb(k_c_offset),(h1b - 503 & 1 + noab * (h5b - 1))) 504 IF (.not.MA_POP_STACK(l_c)) CALL ERRQUIT('ccsdt_o1_2_2',10,MA_ERR) 505 IF (.not.MA_POP_STACK(l_c_sort)) CALL ERRQUIT('ccsdt_o1_2_2',11,MA 506 &_ERR) 507 END IF 508 END IF 509 END IF 510 next = NXTASK(nprocs,1) 511 END IF 512 count = count + 1 513 END DO 514 END DO 515 next = NXTASK(-nprocs,1) 516 call GA_SYNC() 517 RETURN 518 END 519 SUBROUTINE ccsdt_o1_3(d_a,k_a_offset,d_b,k_b_offset,d_c,k_c_offset 520 &) 521C $Id$ 522C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 523C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 524C i0 ( p2 h1 )_to + = 1 * Sum ( p3 ) * o ( p2 p3 )_o * t ( p3 h1 )_t 525 IMPLICIT NONE 526#include "global.fh" 527#include "mafdecls.fh" 528#include "sym.fh" 529#include "errquit.fh" 530#include "tce.fh" 531 INTEGER d_a 532 INTEGER k_a_offset 533 INTEGER d_b 534 INTEGER k_b_offset 535 INTEGER d_c 536 INTEGER k_c_offset 537 INTEGER NXTASK 538 INTEGER next 539 INTEGER nprocs 540 INTEGER count 541 INTEGER p2b 542 INTEGER h1b 543 INTEGER dimc 544 INTEGER l_c_sort 545 INTEGER k_c_sort 546 INTEGER p3b 547 INTEGER p2b_1 548 INTEGER p3b_1 549 INTEGER p3b_2 550 INTEGER h1b_2 551 INTEGER dim_common 552 INTEGER dima_sort 553 INTEGER dima 554 INTEGER dimb_sort 555 INTEGER dimb 556 INTEGER l_a_sort 557 INTEGER k_a_sort 558 INTEGER l_a 559 INTEGER k_a 560 INTEGER l_b_sort 561 INTEGER k_b_sort 562 INTEGER l_b 563 INTEGER k_b 564 INTEGER l_c 565 INTEGER k_c 566 EXTERNAL NXTASK 567 nprocs = GA_NNODES() 568 count = 0 569 next = NXTASK(nprocs,1) 570 DO p2b = noab+1,noab+nvab 571 DO h1b = 1,noab 572 IF (next.eq.count) THEN 573 IF ((.not.restricted).or.(int_mb(k_spin+p2b-1)+int_mb(k_spin+h1b-1 574 &).ne.4)) THEN 575 IF (int_mb(k_spin+p2b-1) .eq. int_mb(k_spin+h1b-1)) THEN 576 IF (ieor(int_mb(k_sym+p2b-1),int_mb(k_sym+h1b-1)) .eq. ieor(irrep_ 577 &t,irrep_o)) THEN 578 dimc = int_mb(k_range+p2b-1) * int_mb(k_range+h1b-1) 579 IF (.not.MA_PUSH_GET(mt_dbl,dimc,'noname',l_c_sort,k_c_sort)) CALL 580 & ERRQUIT('ccsdt_o1_3',0,MA_ERR) 581 CALL DFILL(dimc,0.0d0,dbl_mb(k_c_sort),1) 582 DO p3b = noab+1,noab+nvab 583 IF (int_mb(k_spin+p2b-1) .eq. int_mb(k_spin+p3b-1)) THEN 584 IF (ieor(int_mb(k_sym+p2b-1),int_mb(k_sym+p3b-1)) .eq. irrep_o) TH 585 &EN 586 CALL TCE_RESTRICTED_2(p2b,p3b,p2b_1,p3b_1) 587 CALL TCE_RESTRICTED_2(p3b,h1b,p3b_2,h1b_2) 588 dim_common = int_mb(k_range+p3b-1) 589 dima_sort = int_mb(k_range+p2b-1) 590 dima = dim_common * dima_sort 591 dimb_sort = int_mb(k_range+h1b-1) 592 dimb = dim_common * dimb_sort 593 IF ((dima .gt. 0) .and. (dimb .gt. 0)) THEN 594 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a_sort,k_a_sort)) CALL 595 & ERRQUIT('ccsdt_o1_3',1,MA_ERR) 596 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a,k_a)) CALL ERRQUIT(' 597 &ccsdt_o1_3',2,MA_ERR) 598 CALL GET_HASH_BLOCK(d_a,dbl_mb(k_a),dima,int_mb(k_a_offset),(p3b_1 599 & - 1 + (noab+nvab) * (p2b_1 - 1))) 600 CALL TCE_SORT_2(dbl_mb(k_a),dbl_mb(k_a_sort),int_mb(k_range+p2b-1) 601 &,int_mb(k_range+p3b-1),1,2,1.0d0) 602 IF (.not.MA_POP_STACK(l_a)) CALL ERRQUIT('ccsdt_o1_3',3,MA_ERR) 603 IF (.not.MA_PUSH_GET(mt_dbl,dimb,'noname',l_b_sort,k_b_sort)) CALL 604 & ERRQUIT('ccsdt_o1_3',4,MA_ERR) 605 IF (.not.MA_PUSH_GET(mt_dbl,dimb,'noname',l_b,k_b)) CALL ERRQUIT(' 606 &ccsdt_o1_3',5,MA_ERR) 607 CALL GET_HASH_BLOCK(d_b,dbl_mb(k_b),dimb,int_mb(k_b_offset),(h1b_2 608 & - 1 + noab * (p3b_2 - noab - 1))) 609 CALL TCE_SORT_2(dbl_mb(k_b),dbl_mb(k_b_sort),int_mb(k_range+p3b-1) 610 &,int_mb(k_range+h1b-1),2,1,1.0d0) 611 IF (.not.MA_POP_STACK(l_b)) CALL ERRQUIT('ccsdt_o1_3',6,MA_ERR) 612 CALL DGEMM('T','N',dima_sort,dimb_sort,dim_common,1.0d0,dbl_mb(k_a 613 &_sort),dim_common,dbl_mb(k_b_sort),dim_common,1.0d0,dbl_mb(k_c_sor 614 &t),dima_sort) 615 IF (.not.MA_POP_STACK(l_b_sort)) CALL ERRQUIT('ccsdt_o1_3',7,MA_ER 616 &R) 617 IF (.not.MA_POP_STACK(l_a_sort)) CALL ERRQUIT('ccsdt_o1_3',8,MA_ER 618 &R) 619 END IF 620 END IF 621 END IF 622 END DO 623 IF (.not.MA_PUSH_GET(mt_dbl,dimc,'noname',l_c,k_c)) CALL ERRQUIT(' 624 &ccsdt_o1_3',9,MA_ERR) 625 CALL TCE_SORT_2(dbl_mb(k_c_sort),dbl_mb(k_c),int_mb(k_range+h1b-1) 626 &,int_mb(k_range+p2b-1),2,1,1.0d0) 627 CALL ADD_HASH_BLOCK(d_c,dbl_mb(k_c),dimc,int_mb(k_c_offset),(h1b - 628 & 1 + noab * (p2b - noab - 1))) 629 IF (.not.MA_POP_STACK(l_c)) CALL ERRQUIT('ccsdt_o1_3',10,MA_ERR) 630 IF (.not.MA_POP_STACK(l_c_sort)) CALL ERRQUIT('ccsdt_o1_3',11,MA_E 631 &RR) 632 END IF 633 END IF 634 END IF 635 next = NXTASK(nprocs,1) 636 END IF 637 count = count + 1 638 END DO 639 END DO 640 next = NXTASK(-nprocs,1) 641 call GA_SYNC() 642 RETURN 643 END 644 SUBROUTINE ccsdt_o1_4(d_a,k_a_offset,d_b,k_b_offset,d_c,k_c_offset 645 &) 646C $Id$ 647C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 648C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 649C i0 ( p2 h1 )_to + = 1 * Sum ( p4 h3 ) * o ( h3 p4 )_o * t ( p2 p4 h1 h3 )_t 650 IMPLICIT NONE 651#include "global.fh" 652#include "mafdecls.fh" 653#include "sym.fh" 654#include "errquit.fh" 655#include "tce.fh" 656 INTEGER d_a 657 INTEGER k_a_offset 658 INTEGER d_b 659 INTEGER k_b_offset 660 INTEGER d_c 661 INTEGER k_c_offset 662 INTEGER NXTASK 663 INTEGER next 664 INTEGER nprocs 665 INTEGER count 666 INTEGER p2b 667 INTEGER h1b 668 INTEGER dimc 669 INTEGER l_c_sort 670 INTEGER k_c_sort 671 INTEGER h3b 672 INTEGER p4b 673 INTEGER h3b_1 674 INTEGER p4b_1 675 INTEGER p2b_2 676 INTEGER p4b_2 677 INTEGER h1b_2 678 INTEGER h3b_2 679 INTEGER dim_common 680 INTEGER dima_sort 681 INTEGER dima 682 INTEGER dimb_sort 683 INTEGER dimb 684 INTEGER l_a_sort 685 INTEGER k_a_sort 686 INTEGER l_a 687 INTEGER k_a 688 INTEGER l_b_sort 689 INTEGER k_b_sort 690 INTEGER l_b 691 INTEGER k_b 692 INTEGER l_c 693 INTEGER k_c 694 EXTERNAL NXTASK 695 nprocs = GA_NNODES() 696 count = 0 697 next = NXTASK(nprocs,1) 698 DO p2b = noab+1,noab+nvab 699 DO h1b = 1,noab 700 IF (next.eq.count) THEN 701 IF ((.not.restricted).or.(int_mb(k_spin+p2b-1)+int_mb(k_spin+h1b-1 702 &).ne.4)) THEN 703 IF (int_mb(k_spin+p2b-1) .eq. int_mb(k_spin+h1b-1)) THEN 704 IF (ieor(int_mb(k_sym+p2b-1),int_mb(k_sym+h1b-1)) .eq. ieor(irrep_ 705 &t,irrep_o)) THEN 706 dimc = int_mb(k_range+p2b-1) * int_mb(k_range+h1b-1) 707 IF (.not.MA_PUSH_GET(mt_dbl,dimc,'noname',l_c_sort,k_c_sort)) CALL 708 & ERRQUIT('ccsdt_o1_4',0,MA_ERR) 709 CALL DFILL(dimc,0.0d0,dbl_mb(k_c_sort),1) 710 DO h3b = 1,noab 711 DO p4b = noab+1,noab+nvab 712 IF (int_mb(k_spin+h3b-1) .eq. int_mb(k_spin+p4b-1)) THEN 713 IF (ieor(int_mb(k_sym+h3b-1),int_mb(k_sym+p4b-1)) .eq. irrep_o) TH 714 &EN 715 CALL TCE_RESTRICTED_2(h3b,p4b,h3b_1,p4b_1) 716 CALL TCE_RESTRICTED_4(p2b,p4b,h1b,h3b,p2b_2,p4b_2,h1b_2,h3b_2) 717 dim_common = int_mb(k_range+h3b-1) * int_mb(k_range+p4b-1) 718 dima_sort = 1 719 dima = dim_common * dima_sort 720 dimb_sort = int_mb(k_range+p2b-1) * int_mb(k_range+h1b-1) 721 dimb = dim_common * dimb_sort 722 IF ((dima .gt. 0) .and. (dimb .gt. 0)) THEN 723 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a_sort,k_a_sort)) CALL 724 & ERRQUIT('ccsdt_o1_4',1,MA_ERR) 725 IF (.not.MA_PUSH_GET(mt_dbl,dima,'noname',l_a,k_a)) CALL ERRQUIT(' 726 &ccsdt_o1_4',2,MA_ERR) 727 CALL GET_HASH_BLOCK(d_a,dbl_mb(k_a),dima,int_mb(k_a_offset),(p4b_1 728 & - 1 + (noab+nvab) * (h3b_1 - 1))) 729 CALL TCE_SORT_2(dbl_mb(k_a),dbl_mb(k_a_sort),int_mb(k_range+h3b-1) 730 &,int_mb(k_range+p4b-1),2,1,1.0d0) 731 IF (.not.MA_POP_STACK(l_a)) CALL ERRQUIT('ccsdt_o1_4',3,MA_ERR) 732 IF (.not.MA_PUSH_GET(mt_dbl,dimb,'noname',l_b_sort,k_b_sort)) CALL 733 & ERRQUIT('ccsdt_o1_4',4,MA_ERR) 734 IF (.not.MA_PUSH_GET(mt_dbl,dimb,'noname',l_b,k_b)) CALL ERRQUIT(' 735 &ccsdt_o1_4',5,MA_ERR) 736 IF ((p4b .lt. p2b) .and. (h3b .lt. h1b)) THEN 737 CALL GET_HASH_BLOCK(d_b,dbl_mb(k_b),dimb,int_mb(k_b_offset),(h1b_2 738 & - 1 + noab * (h3b_2 - 1 + noab * (p2b_2 - noab - 1 + nvab * (p4b_ 739 &2 - noab - 1))))) 740 CALL TCE_SORT_4(dbl_mb(k_b),dbl_mb(k_b_sort),int_mb(k_range+p4b-1) 741 &,int_mb(k_range+p2b-1),int_mb(k_range+h3b-1),int_mb(k_range+h1b-1) 742 &,4,2,1,3,1.0d0) 743 END IF 744 IF ((p4b .lt. p2b) .and. (h1b .le. h3b)) THEN 745 CALL GET_HASH_BLOCK(d_b,dbl_mb(k_b),dimb,int_mb(k_b_offset),(h3b_2 746 & - 1 + noab * (h1b_2 - 1 + noab * (p2b_2 - noab - 1 + nvab * (p4b_ 747 &2 - noab - 1))))) 748 CALL TCE_SORT_4(dbl_mb(k_b),dbl_mb(k_b_sort),int_mb(k_range+p4b-1) 749 &,int_mb(k_range+p2b-1),int_mb(k_range+h1b-1),int_mb(k_range+h3b-1) 750 &,3,2,1,4,-1.0d0) 751 END IF 752 IF ((p2b .le. p4b) .and. (h3b .lt. h1b)) THEN 753 CALL GET_HASH_BLOCK(d_b,dbl_mb(k_b),dimb,int_mb(k_b_offset),(h1b_2 754 & - 1 + noab * (h3b_2 - 1 + noab * (p4b_2 - noab - 1 + nvab * (p2b_ 755 &2 - noab - 1))))) 756 CALL TCE_SORT_4(dbl_mb(k_b),dbl_mb(k_b_sort),int_mb(k_range+p2b-1) 757 &,int_mb(k_range+p4b-1),int_mb(k_range+h3b-1),int_mb(k_range+h1b-1) 758 &,4,1,2,3,-1.0d0) 759 END IF 760 IF ((p2b .le. p4b) .and. (h1b .le. h3b)) THEN 761 CALL GET_HASH_BLOCK(d_b,dbl_mb(k_b),dimb,int_mb(k_b_offset),(h3b_2 762 & - 1 + noab * (h1b_2 - 1 + noab * (p4b_2 - noab - 1 + nvab * (p2b_ 763 &2 - noab - 1))))) 764 CALL TCE_SORT_4(dbl_mb(k_b),dbl_mb(k_b_sort),int_mb(k_range+p2b-1) 765 &,int_mb(k_range+p4b-1),int_mb(k_range+h1b-1),int_mb(k_range+h3b-1) 766 &,3,1,2,4,1.0d0) 767 END IF 768 IF (.not.MA_POP_STACK(l_b)) CALL ERRQUIT('ccsdt_o1_4',6,MA_ERR) 769 CALL DGEMM('T','N',dima_sort,dimb_sort,dim_common,1.0d0,dbl_mb(k_a 770 &_sort),dim_common,dbl_mb(k_b_sort),dim_common,1.0d0,dbl_mb(k_c_sor 771 &t),dima_sort) 772 IF (.not.MA_POP_STACK(l_b_sort)) CALL ERRQUIT('ccsdt_o1_4',7,MA_ER 773 &R) 774 IF (.not.MA_POP_STACK(l_a_sort)) CALL ERRQUIT('ccsdt_o1_4',8,MA_ER 775 &R) 776 END IF 777 END IF 778 END IF 779 END DO 780 END DO 781 IF (.not.MA_PUSH_GET(mt_dbl,dimc,'noname',l_c,k_c)) CALL ERRQUIT(' 782 &ccsdt_o1_4',9,MA_ERR) 783 CALL TCE_SORT_2(dbl_mb(k_c_sort),dbl_mb(k_c),int_mb(k_range+h1b-1) 784 &,int_mb(k_range+p2b-1),2,1,1.0d0) 785 CALL ADD_HASH_BLOCK(d_c,dbl_mb(k_c),dimc,int_mb(k_c_offset),(h1b - 786 & 1 + noab * (p2b - noab - 1))) 787 IF (.not.MA_POP_STACK(l_c)) CALL ERRQUIT('ccsdt_o1_4',10,MA_ERR) 788 IF (.not.MA_POP_STACK(l_c_sort)) CALL ERRQUIT('ccsdt_o1_4',11,MA_E 789 &RR) 790 END IF 791 END IF 792 END IF 793 next = NXTASK(nprocs,1) 794 END IF 795 count = count + 1 796 END DO 797 END DO 798 next = NXTASK(-nprocs,1) 799 call GA_SYNC() 800 RETURN 801 END 802