1c 2c ============================================== 3c Create effective Hamiltonian 4c ============================================== 5c 6 subroutine tce_heff(d_em,k_e_offsetm,k_r1_offsetm, 7 1 k_r2_offsetm,k_r3_offsetm,k_r4_offsetm,d_r1m,d_r2m,d_r3m,d_r4m, 8 2 needt1,needt2,needr3act,needr4act,rtdb) 9 implicit none 10#include "tce.fh" 11#include "mafdecls.fh" 12#include "stdio.fh" 13#include "rtdb.fh" 14#include "errquit.fh" 15#include "sym.fh" 16#include "tce_mrcc.fh" 17#include "global.fh" 18#include "tce_main.fh" 19 20 integer d_em(maxref) 21 integer k_e_offsetm(maxref) 22 integer iref,i 23 double precision corr 24 logical nodezero 25 integer k_r1_offsetm(maxref) 26 integer k_r2_offsetm(maxref) 27 integer k_r3_offsetm(maxref) 28 integer k_r4_offsetm(maxref) 29 integer d_r1m(maxref),d_r2m(maxref) 30 integer d_r3m(maxref),d_r4m(maxref) 31 logical needt1,needt2,needr3act 32 logical needr4act 33 integer rtdb 34 35 nodezero = (ga_nodeid().eq.0) 36 37 if(lusesub)call ga_zero(g_heff) 38 39 do i=1,nref*nref 40 dbl_mb(k_heff+i-1) = 0.0d0 41 enddo 42 43 do iref=1,nref 44 45 if((int_mb(k_refafi+iref-1).eq.int_mb(k_innodes 46 1 +ga_nnodes()+ga_nodeid())).or.(.not.lusesub)) then 47 48 call get_block(d_em(iref),corr,1,0) 49 dbl_mb(k_heff+iref-1+(iref-1)*nref) = corr+duhfens(iref) 50 51 if(lusesub) then 52 53c write(6,*)ga_nodeid(),corr+duhfens(iref),iref 54 call ga_put(g_heff,nref*(iref-1)+iref,nref*(iref-1)+iref,1,1, 55 1 corr+duhfens(iref),1) 56 endif 57 58 endif 59 60 enddo 61 62 call tce_heff_offdiagonal(k_r1_offsetm,k_r2_offsetm, 63 1 k_r3_offsetm,k_r4_offsetm,d_r1m,d_r2m,d_r3m,d_r4m,needt1, 64 2 needt2,needr3act,needr4act,rtdb) 65 66c if(nodezero) then 67c call ma_print(dbl_mb(k_heff),nref,nref,'Heff') 68c endif 69 70 return 71 end 72c 73c ============================================== 74c Add offdiagonal elements of Heff 75c ============================================== 76c 77 subroutine tce_heff_offdiagonal(k_r1_offsetm, 78 1 k_r2_offsetm,k_r3_offsetm,k_r4_offsetm,d_r1m,d_r2m,d_r3m, 79 2 d_r4m,needt1,needt2,needr3act,needr4act,rtdb) 80 implicit none 81#include "tce.fh" 82#include "mafdecls.fh" 83#include "stdio.fh" 84#include "rtdb.fh" 85#include "errquit.fh" 86#include "sym.fh" 87#include "tce_mrcc.fh" 88#include "global.fh" 89#include "tce_main.fh" 90 91 integer rtdb 92 logical nodezero 93 integer k_r1_offsetm(maxref) 94 integer k_r2_offsetm(maxref) 95 integer k_r3_offsetm(maxref) 96 integer k_r4_offsetm(maxref) 97 integer d_r1m(maxref),d_r2m(maxref) 98 integer d_r3m(maxref) 99 integer d_r4m(maxref) 100 integer iref 101 integer i,j,p5b,h6b,k 102 integer size,l,m,n,o 103 integer l_r1,k_r1,l_r2,k_r2 104 integer l_r3,k_r3 105 integer l_r4,k_r4 106 integer p1b,p2b,h3b,h4b 107 integer p3b,p7b,p8b 108 INTEGER p4b 109 INTEGER p6b 110 INTEGER h1b 111 INTEGER h2b 112 integer h1,h2,h3,p4,p5,p6 113 integer h4,p7,p8 114 integer orbindex(8) 115 integer orbspin(8) 116 integer ioccnew(maxorb,2) 117 integer iu 118 logical isfound 119 logical needt1,needt2,needr3act 120 logical needr4act 121c logical lusescffv 122 integer is,iocc0(maxorb,2) 123 integer i1,i2,dist 124 double precision dsmult 125c logical limprovet 126 EXTERNAL NXTASKsub 127 EXTERNAL NXTASK 128 INTEGER NXTASKsub 129 INTEGER NXTASK 130 INTEGER nxt 131 INTEGER nprocs 132 INTEGER count 133 134 nodezero = (ga_nodeid().eq.0) 135 isfound = .false. 136 137c if (.not.rtdb_get(rtdb,'mrcc:usescffermiv',mt_log,1,lusescffv)) 138c 1 lusescffv = .false. 139c if (.not.rtdb_get(rtdb,'mrcc:improvetiling',mt_log,1,limprovet)) 140c 1 limprovet = .false. 141 142 do is=1,2 143 do i=1,nmo(is) 144 if(i.le.nocc(is)) then 145 iocc0(i,is) = 1 146 else 147 iocc0(i,is) = 0 148 end if 149 end do 150 end do 151 152 do iref=1,nref 153 154 if((int_mb(k_refafi+iref-1).eq.int_mb(k_innodes 155 1 +ga_nnodes()+ga_nodeid())).or.(.not.lusesub)) then 156 157 k_sym = k_symm(iref) 158 k_offset = k_offsetm(iref) 159 k_range = k_rangem(iref) 160 k_spin = k_spinm(iref) 161 k_movecs_sorted = k_movecs_sortedm(iref) 162 k_active = k_active_tmpm(iref) 163 164 noa = nblcks(1,iref) 165 nob = nblcks(2,iref) 166 nva = nblcks(3,iref) 167 nvb = nblcks(4,iref) 168 169 noab = noa+nob 170 nvab = nva+nvb 171c 172c --------------- 173c R1 active 174c --------------- 175c 176 do p5b = noab+1,noab+nvab 177 do h6b = 1,noab 178 179 k = 0 180 181 if (int_mb(k_spin+p5b-1) .eq. int_mb(k_spin+h6b-1)) then 182 if (ieor(int_mb(k_sym+p5b-1),int_mb(k_sym+h6b-1)).eq.irrep_t)then 183 if ((.not.restricted).or.(int_mb(k_spin+p5b-1)+int_mb(k_spin+h6b-1 184 &).ne.4)) then 185 if(log_mb(k_isactive(iref)+p5b-1).and. 186 &log_mb(k_isactive(iref)+h6b-1)) then 187 188 size = int_mb(k_range+h6b-1) * int_mb(k_range+p5b-1) 189 190 if (.not.ma_push_get(mt_dbl,size,'r1mi',l_r1,k_r1)) 191 1 call errquit('tce_mrcc_iface_r1: MA problem',0,MA_ERR) 192 193 call get_hash_block(d_r1m(iref),dbl_mb(k_r1),size, 194 1 int_mb(k_r1_offsetm(iref)),h6b-1+noab*(p5b-noab-1)) 195 196 k=0 197 do i=1,int_mb(k_range+p5b-1) 198 do j=1,int_mb(k_range+h6b-1) 199 k = k + 1 200 201 orbspin(1) = int_mb(k_spin+p5b-1)-1 202 orbspin(2) = int_mb(k_spin+h6b-1)-1 203 204 orbindex(1) = (1 - orbspin(1)+ 205 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p5b-1)+i-1))/2 206 orbindex(2) = (1 - orbspin(2)+ 207 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h6b-1)+j-1))/2 208 209 orbindex(1) = moindexes(orbindex(1),orbspin(1)+1,iref) 210 orbindex(2) = moindexes(orbindex(2),orbspin(2)+1,iref) 211 212cjb =========================================================================== 213 214 if(isactive(orbindex(1),orbspin(1)+1).and. 215 1 isactive(orbindex(2),orbspin(2)+1).or.(.not.limprovet)) then 216 217ccc ===== 218c if(nodezero)write(6,"('ACTIVITY: ',2L2)") 219c 1 isactive(orbindex(1),orbspin(1)+1), 220c 1 isactive(orbindex(2),orbspin(2)+1) 221c if(nodezero)write(6,"('DEBUG: ',5I4)") 222c 1 orbindex(1),orbspin(1), 223c 1 orbindex(2),orbspin(2),iref 224ccc ===== 225 226 do iu=1,2 227 do n=1,nmo(iu) 228 ioccnew(n,iu) = iocc(n,iref,iu) 229 enddo 230 enddo 231 232 if(iocc(orbindex(1),iref,orbspin(1)+1).eq. 233 1 iocc(orbindex(2),iref,orbspin(2)+1)) then 234 goto 200 235 endif 236 237 ioccnew(orbindex(1),orbspin(1)+1) = iocc(orbindex(2),iref, 238 1 orbspin(2)+1) 239 ioccnew(orbindex(2),orbspin(2)+1) = iocc(orbindex(1),iref, 240 1 orbspin(1)+1) 241 242 do n=1,nref 243 isfound = .true. 244 do iu=1,2 245 do o=1,nmo(iu) 246 isfound = isfound.and.(iocc(o,n,iu).eq.ioccnew(o,iu)) 247 enddo 248 enddo 249 if(isfound) then 250c write(LuOut,"('Internal amplitude',I4,'->',I4)")iref,n 251c write(LuOut,"('1Internal amplitude',I4,'->',I4,2F16.12)")iref,n, 252c 1 dbl_mb(k_r1+m-1) 253 if(iref.ne.n) then 254 255 dist = 0 256 do iu=1,2 257 do i1=1,nmo(iu) 258 ioccnew(i1,iu) = iocc(i1,iref,iu) 259 enddo 260 enddo 261 262 i2 = 0 263 do i1=min(orbindex(1),orbindex(2)), 264 1 max(orbindex(1),orbindex(2)) 265 i2 = i2 + 1 266 if(i2.lt.abs(orbindex(1)-orbindex(2))) then 267 if(iocc(i1+1,iref,orbspin(1)+1).eq.1) then 268 dist=dist+1 269 endif 270 endif 271 enddo 272 273 dsmult = 1.0d0 274 275 if(mod(dist,2).eq.1) dsmult = -dsmult 276 277 dbl_mb(k_heff+n-1+(iref-1)*nref) = dbl_mb(k_r1+k-1)*dsmult 278c 1 dbl_mb(k_heff+n-1+(iref-1)*nref) 279 280 if(lusesub) then 281 call ga_put(g_heff,nref*(iref-1)+n,nref*(iref-1)+n,1,1, 282 1 dbl_mb(k_r1+k-1)*dsmult,1) 283 endif 284 285 endif 286 goto 200 287 endif 288 enddo 289 290 200 continue 291 292 if((.not.isfound).and.(abs(dbl_mb(k_r1+k-1)).gt.1.0d-5)) then 293 if(nodezero) then 294 write(LuOut,"('DEBUG: ',4F16.12)")dbl_mb(k_r1+k-1) 295 write(LuOut,"('YOU ARE USING INCOMPLETE MODEL SPACE!')") 296 endif 297c call errquit('YOU ARE USING INCOMPLETE MODEL SPACE!',1,MA_ERR) 298 endif 299 300 endif ! active 301 302 enddo 303 enddo 304 305 if (.not.ma_pop_stack(l_r1)) 306 1 call errquit('tce_mrcc_iface_r1: MA problem',1,MA_ERR) 307 308 endif 309 endif 310 endif 311 endif 312 enddo 313 enddo 314c 315c --------------- 316c R2 active 317c --------------- 318c 319 nxt = 0 320 count = 0 321 if(limprovet) then 322 if(lusesub) then 323 nprocs=GA_pgroup_NNODES(mypgid) 324 nxt=NXTASKsub(nprocs,1,mypgid) 325 else 326 nprocs=GA_NNODES() 327 nxt=NXTASK(nprocs,1) 328 endif 329 count = 0 330 endif 331 332 DO p1b = noab+1,noab+nvab 333 DO p2b = p1b,noab+nvab 334 DO h3b = 1,noab 335 DO h4b = h3b,noab 336 337 IF ((nxt.eq.count).or.(.not.limprovet)) THEN 338 339 IF (int_mb(k_spin+p1b-1)+int_mb(k_spin+p2b-1) .eq. int_mb(k_spin+h 340 &3b-1)+int_mb(k_spin+h4b-1)) THEN 341 342 IF (ieor(int_mb(k_sym+p1b-1),ieor(int_mb(k_sym+p2b-1),ieor(int_mb( 343 &k_sym+h3b-1),int_mb(k_sym+h4b-1)))) .eq. irrep_t) THEN 344 345 IF ((.not.restricted).or.(int_mb(k_spin+p1b-1)+int_mb(k_spin+p2b-1 346 &)+int_mb(k_spin+h3b-1)+int_mb(k_spin+h4b-1).ne.8)) THEN 347 348 if(log_mb(k_isactive(iref)+p1b-1).and. 349 1 log_mb(k_isactive(iref)+p2b-1).and. 350 2 log_mb(k_isactive(iref)+h3b-1).and. 351 3 log_mb(k_isactive(iref)+h4b-1)) then 352 353 size = int_mb(k_range+p1b-1) * int_mb(k_range+p2b-1) * int_ 354 &mb(k_range+h3b-1) * int_mb(k_range+h4b-1) 355 356 if (.not.ma_push_get(mt_dbl,size,'r2mi',l_r2,k_r2)) 357 1 call errquit('tce_mrcc_iface_r2: MA problem',0,MA_ERR) 358 359 call get_hash_block(d_r2m(iref),dbl_mb(k_r2),size, 360 1 int_mb(k_r2_offsetm(iref)),h4b-1+noab*(h3b-1+noab*(p2b- 361 &noab-1+nvab*(p1b - noab - 1)))) 362 363c write(LuOut,"(I4,L3,L3,L3,L3)") 364c 1 iref,log_mb(k_isactive(iref)+p1b-1), 365c 1 log_mb(k_isactive(iref)+p2b-1),log_mb(k_isactive(iref)+h3b-1), 366c 1 log_mb(k_isactive(iref)+h4b-1) 367 m = 0 368 369 do i=1,int_mb(k_range+p1b-1) 370 do j=1,int_mb(k_range+p2b-1) 371 do k=1,int_mb(k_range+h3b-1) 372 do l=1,int_mb(k_range+h4b-1) 373 m = m + 1 374c write(LuOut,"(I4,'(',I4,I4,I4,I4,'):',2F16.12)") 375c 1 iref,i,j,k,l,dbl_mb(k_r2+m-1) 376c write(LuOut,*)int_mb(k_spin+p1b-1) 377 378 orbspin(1) = int_mb(k_spin+p1b-1)-1 379 orbspin(2) = int_mb(k_spin+p2b-1)-1 380 orbspin(3) = int_mb(k_spin+h3b-1)-1 381 orbspin(4) = int_mb(k_spin+h4b-1)-1 382 383 orbindex(1) = (1 - orbspin(1)+ 384 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p1b-1)+i-1))/2 385 orbindex(2) = (1 - orbspin(2)+ 386 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p2b-1)+j-1))/2 387 orbindex(3) = (1 - orbspin(3)+ 388 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h3b-1)+k-1))/2 389 orbindex(4) = (1 - orbspin(4)+ 390 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h4b-1)+l-1))/2 391 392 orbindex(1) = moindexes(orbindex(1),orbspin(1)+1,iref) 393 orbindex(2) = moindexes(orbindex(2),orbspin(2)+1,iref) 394 orbindex(3) = moindexes(orbindex(3),orbspin(3)+1,iref) 395 orbindex(4) = moindexes(orbindex(4),orbspin(4)+1,iref) 396 397cjb =========================================================================== 398 399 if(isactive(orbindex(1),orbspin(1)+1).and. 400 1 isactive(orbindex(2),orbspin(2)+1).and. 401 2 isactive(orbindex(3),orbspin(3)+1).and. 402 3 isactive(orbindex(4),orbspin(4)+1).or.(.not.limprovet)) then 403 404c write(LuOut,"('Real indexes: [',I4,I4,I4,I4,']', 405c 1 '[',I2,I2,I2,I2,']')") 406c 1 orbindex(1),orbindex(2),orbindex(3),orbindex(4), 407c 1 orbspin(1),orbspin(2),orbspin(3),orbspin(4) 408 409 do iu=1,2 410 do n=1,nmo(iu) 411 ioccnew(n,iu) = iocc(n,iref,iu) 412 enddo 413 enddo 414 415 if(((orbindex(1).eq.orbindex(2)).and.(orbspin(1).eq.orbspin(2))) 416 1 .or.((orbindex(3).eq.orbindex(4)).and.(orbspin(3).eq. 417 2 orbspin(4)))) then 418 goto 100 419 endif 420 421 ioccnew(orbindex(1),orbspin(1)+1) = iocc(orbindex(3),iref, 422 1 orbspin(3)+1) 423 ioccnew(orbindex(2),orbspin(2)+1) = iocc(orbindex(4),iref, 424 1 orbspin(4)+1) 425 ioccnew(orbindex(3),orbspin(3)+1) = iocc(orbindex(1),iref, 426 1 orbspin(1)+1) 427 ioccnew(orbindex(4),orbspin(4)+1) = iocc(orbindex(2),iref, 428 1 orbspin(2)+1) 429 430 do n=1,nref 431 isfound = .true. 432 do iu=1,2 433 do o=1,nmo(iu) 434 isfound = isfound.and.(iocc(o,n,iu).eq.ioccnew(o,iu)) 435 enddo 436 enddo 437 if(isfound) then 438c write(LuOut,"('2Internal amplitude',I4,'->',I4,2F16.12)")iref,n, 439c 1 dbl_mb(k_r2+m-1) 440 if(iref.ne.n) then 441 442 dist = 0 443 do iu=1,2 444 do i1=1,nmo(iu) 445 ioccnew(i1,iu) = iocc(i1,iref,iu) 446 enddo 447 enddo 448 449 i2 = 0 450 do i1=min(orbindex(1),orbindex(3)), 451 1 max(orbindex(1),orbindex(3)) 452 i2 = i2 + 1 453 if(i2.lt.abs(orbindex(1)-orbindex(3))) then 454 if(iocc(i1+1,iref,orbspin(1)+1).eq.1) then 455 dist=dist+1 456 endif 457 endif 458 enddo 459 460 ioccnew(orbindex(1),orbspin(1)+1) = iocc(orbindex(3),iref, 461 1 orbspin(3)+1) 462 ioccnew(orbindex(3),orbspin(3)+1) = iocc(orbindex(1),iref, 463 1 orbspin(1)+1) 464 465 i2 = 0 466 do i1=min(orbindex(2),orbindex(4)), 467 1 max(orbindex(2),orbindex(4)) 468 i2 = i2 + 1 469 if(i2.lt.abs(orbindex(2)-orbindex(4))) then 470 if(ioccnew(i1+1,orbspin(2)+1).eq.1) then 471 dist=dist+1 472 endif 473 endif 474 enddo 475 476 dsmult = 1.0d0 477 if(mod(dist,2).eq.1) dsmult = -dsmult 478 479 480 dbl_mb(k_heff+n-1+(iref-1)*nref) = dbl_mb(k_r2+m-1)*dsmult 481c 1 dbl_mb(k_heff+n-1+(iref-1)*nref) 482 483 if(lusesub) then 484 call ga_put(g_heff,nref*(iref-1)+n,nref*(iref-1)+n,1,1, 485 1 dbl_mb(k_r2+m-1)*dsmult,1) 486 endif 487 488 489 endif 490 goto 100 491 endif 492 enddo 493 494 100 continue 495 496 if((.not.isfound).and.(abs(dbl_mb(k_r2+m-1)).gt.1.0d-5)) then 497 if(nodezero) then 498 write(LuOut,"('DEBUG: ',4F16.12)")dbl_mb(k_r2+m-1) 499 write(LuOut,"('YOU ARE USING INCOMPLETE MODEL SPACE!')") 500 endif 501c call errquit('YOU ARE USING INCOMPLETE MODEL SPACE!',2,MA_ERR) 502 endif 503 504 endif ! act 505 506 enddo 507 enddo 508 enddo 509 enddo 510 511 if (.not.ma_pop_stack(l_r2)) 512 1 call errquit('tce_mrcc_iface_r2: MA problem',1,MA_ERR) 513 END IF 514 END IF 515 END IF 516 endif 517 if(limprovet) then 518 if(lusesub) then 519 nxt=NXTASKsub(nprocs,1,mypgid) 520 else 521 nxt=NXTASK(nprocs,1) 522 endif 523 endif 524 endif 525 if(limprovet)count = count + 1 526 END DO 527 END DO 528 END DO 529 END DO 530 531 if(limprovet) then 532 if(lusesub) then 533 nxt=NXTASKsub(-nprocs,1,mypgid) 534 call GA_Pgroup_SYNC(mypgid) 535 else 536 nxt=NXTASK(-nprocs,1) 537 call GA_SYNC() 538 endif 539 endif 540c 541c --------------- 542c R3 active 543c --------------- 544c 545 if(needr3act) then 546 DO p4b = noab+1,noab+nvab 547 DO p5b = p4b,noab+nvab 548 DO p6b = p5b,noab+nvab 549 DO h1b = 1,noab 550 DO h2b = h1b,noab 551 DO h3b = h2b,noab 552 IF ((.not.restricted).or.(int_mb(k_spin+p4b-1)+int_mb(k_spin+p5b-1 553 &)+int_mb(k_spin+p6b-1)+int_mb(k_spin+h1b-1)+int_mb(k_spin+h2b-1)+i 554 &nt_mb(k_spin+h3b-1).ne.12)) THEN 555 IF (int_mb(k_spin+p4b-1)+int_mb(k_spin+p5b-1)+int_mb(k_spin+p6b-1) 556 & .eq. int_mb(k_spin+h1b-1)+int_mb(k_spin+h2b-1)+int_mb(k_spin+h3b- 557 &1)) THEN 558 IF (ieor(int_mb(k_sym+p4b-1),ieor(int_mb(k_sym+p5b-1),ieor(int_mb( 559 &k_sym+p6b-1),ieor(int_mb(k_sym+h1b-1),ieor(int_mb(k_sym+h2b-1),int 560 &_mb(k_sym+h3b-1)))))) .eq. ieor(irrep_v,irrep_t)) THEN 561 IF ((log_mb(k_active+p4b-1).eqv..true.).and.(log_mb(k_active+p5b-1 562 &).eqv..true.).and.(log_mb(k_active+p6b-1).eqv..true.).and.(log_mb( 563 &k_active+h1b-1).eqv..true.).and.(log_mb(k_active+h2b-1).eqv..true. 564 &).and.(log_mb(k_active+h3b-1).eqv..true.)) THEN 565 566 size = int_mb(k_range+p4b-1) * int_mb(k_range+p5b-1) * int_ 567 &mb(k_range+p6b-1) * int_mb(k_range+h1b-1) * int_mb(k_range+h2b-1) 568 &* int_mb(k_range+h3b-1) 569 570 if (.not.ma_push_get(mt_dbl,size,'r3mi',l_r3,k_r3)) 571 1 call errquit('tce_mrcc_iface_r3: MA problem',0,MA_ERR) 572 573 call get_hash_block(d_r3m(iref),dbl_mb(k_r3),size, 574 1 int_mb(k_r3_offsetm(iref)),h3b - 1 + noab * (h2b - 1 + noab * 575 1 (h1b - 1 + noab * (p6b - noab - 1 + nvab * (p5b - noab - 1 + 576 1 nvab * (p4b - noab - 1)))))) 577 578 m = 0 579 580 orbspin(1) = int_mb(k_spin+p4b-1)-1 581 orbspin(2) = int_mb(k_spin+p5b-1)-1 582 orbspin(3) = int_mb(k_spin+p6b-1)-1 583 orbspin(4) = int_mb(k_spin+h1b-1)-1 584 orbspin(5) = int_mb(k_spin+h2b-1)-1 585 orbspin(6) = int_mb(k_spin+h3b-1)-1 586 587 do p4=1,int_mb(k_range+p4b-1) 588 do p5=1,int_mb(k_range+p5b-1) 589 do p6=1,int_mb(k_range+p6b-1) 590 do h1=1,int_mb(k_range+h1b-1) 591 do h2=1,int_mb(k_range+h2b-1) 592 do h3=1,int_mb(k_range+h3b-1) 593 594 m = m + 1 595 596 orbindex(1) = (1 - orbspin(1)+ 597 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p4b-1)+p4-1))/2 598 orbindex(2) = (1 - orbspin(2)+ 599 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p5b-1)+p5-1))/2 600 orbindex(3) = (1 - orbspin(3)+ 601 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p6b-1)+p6-1))/2 602 orbindex(4) = (1 - orbspin(4)+ 603 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h1b-1)+h1-1))/2 604 orbindex(5) = (1 - orbspin(5)+ 605 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h2b-1)+h2-1))/2 606 orbindex(6) = (1 - orbspin(6)+ 607 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h3b-1)+h3-1))/2 608 609 orbindex(1) = moindexes(orbindex(1),orbspin(1)+1,iref) 610 orbindex(2) = moindexes(orbindex(2),orbspin(2)+1,iref) 611 orbindex(3) = moindexes(orbindex(3),orbspin(3)+1,iref) 612 orbindex(4) = moindexes(orbindex(4),orbspin(4)+1,iref) 613 orbindex(5) = moindexes(orbindex(5),orbspin(5)+1,iref) 614 orbindex(6) = moindexes(orbindex(6),orbspin(6)+1,iref) 615 616c write(LuOut,"('Real indexes: [',I4,I4,I4,I4,I4,I4,']')") 617c 1 orbindex(1),orbindex(2),orbindex(3),orbindex(4),orbindex(5), 618c 1 orbindex(6) 619c write(LuOut,"('Spin indexes : [',I4,I4,I4,I4,I4,I4,']')") 620c 1 orbspin(1),orbspin(2),orbspin(3),orbspin(4),orbspin(5),orbspin(6) 621 622 do iu=1,2 623 do n=1,nbf 624 ioccnew(n,iu) = iocc(n,iref,iu) 625 enddo 626 enddo 627 628 if((iocc(orbindex(1),iref,orbspin(1)+1).eq. 629 1 iocc(orbindex(4),iref,orbspin(4)+1)).or. 630 2 (iocc(orbindex(2),iref,orbspin(2)+1).eq. 631 3 iocc(orbindex(5),iref,orbspin(5)+1)).or. 632 4 (iocc(orbindex(3),iref,orbspin(3)+1).eq. 633 1 iocc(orbindex(6),iref,orbspin(6)+1))) then 634 goto 300 635 endif 636 637 if((orbspin(1).ne.orbspin(4)).or. 638 1 (orbspin(2).ne.orbspin(5)).or. 639 2 (orbspin(3).ne.orbspin(6))) then 640 goto 300 641 endif 642 643 if( 644 1 ((orbindex(1).eq.orbindex(2)).and.(orbspin(1).eq.orbspin(2))).or. 645 1 ((orbindex(1).eq.orbindex(3)).and.(orbspin(1).eq.orbspin(3))).or. 646 1 ((orbindex(2).eq.orbindex(3)).and.(orbspin(2).eq.orbspin(3))).or. 647 1 ((orbindex(4).eq.orbindex(5)).and.(orbspin(4).eq.orbspin(5))).or. 648 1 ((orbindex(4).eq.orbindex(6)).and.(orbspin(4).eq.orbspin(6))).or. 649 1 ((orbindex(5).eq.orbindex(6)).and.(orbspin(5).eq.orbspin(6))) 650 1 ) then 651 goto 300 652 endif 653 654 ioccnew(orbindex(1),orbspin(1)+1) = iocc(orbindex(4),iref, 655 1 orbspin(4)+1) 656 ioccnew(orbindex(2),orbspin(2)+1) = iocc(orbindex(5),iref, 657 1 orbspin(5)+1) 658 ioccnew(orbindex(3),orbspin(3)+1) = iocc(orbindex(6),iref, 659 1 orbspin(6)+1) 660 ioccnew(orbindex(4),orbspin(4)+1) = iocc(orbindex(1),iref, 661 1 orbspin(1)+1) 662 ioccnew(orbindex(5),orbspin(5)+1) = iocc(orbindex(2),iref, 663 1 orbspin(2)+1) 664 ioccnew(orbindex(6),orbspin(6)+1) = iocc(orbindex(3),iref, 665 1 orbspin(3)+1) 666 667 668 do n=1,nref 669 isfound = .true. 670 do iu=1,2 671 do o=1,nbf 672 isfound = isfound.and.(iocc(o,n,iu).eq.ioccnew(o,iu)) 673 enddo 674 enddo 675 if(isfound) then 676c write(LuOut,"('Internal amplitude',I4,'->',I4,2F16.12)")iref,n, 677c 1 dbl_mb(k_r3+m-1) 678 if(iref.ne.n) then 679cJB START 680 dist = 0 681 do iu=1,2 682 do i1=1,nmo(iu) 683 ioccnew(i1,iu) = iocc(i1,iref,iu) 684 enddo 685 enddo 686 687 i2 = 0 688 do i1=min(orbindex(1),orbindex(4)), 689 1 max(orbindex(1),orbindex(4)) 690 i2 = i2 + 1 691 if(i2.lt.abs(orbindex(1)-orbindex(4))) then 692 if(iocc(i1+1,iref,orbspin(1)+1).eq.1) then 693 dist=dist+1 694 endif 695 endif 696 enddo 697 698 ioccnew(orbindex(1),orbspin(1)+1) = iocc(orbindex(4),iref, 699 1 orbspin(4)+1) 700 ioccnew(orbindex(4),orbspin(4)+1) = iocc(orbindex(1),iref, 701 1 orbspin(1)+1) 702 703 i2 = 0 704 do i1=min(orbindex(2),orbindex(5)), 705 1 max(orbindex(2),orbindex(5)) 706 i2 = i2 + 1 707 if(i2.lt.abs(orbindex(2)-orbindex(5))) then 708 if(ioccnew(i1+1,orbspin(5)+1).eq.1) then 709 dist=dist+1 710 endif 711 endif 712 enddo 713 714 ioccnew(orbindex(2),orbspin(2)+1) = iocc(orbindex(5),iref, 715 1 orbspin(5)+1) 716 ioccnew(orbindex(5),orbspin(5)+1) = iocc(orbindex(2),iref, 717 1 orbspin(2)+1) 718 719 i2 = 0 720 do i1=min(orbindex(3),orbindex(6)), 721 1 max(orbindex(3),orbindex(6)) 722 i2 = i2 + 1 723 if(i2.lt.abs(orbindex(3)-orbindex(6))) then 724 if(ioccnew(i1+1,orbspin(6)+1).eq.1) then 725 dist=dist+1 726 endif 727 endif 728 enddo 729 730 dsmult = 1.0d0 731 if(mod(dist,2).eq.1) dsmult = -dsmult 732 733cJB END 734 dbl_mb(k_heff+n-1+(iref-1)*nref) = 735 1 dbl_mb(k_r3+m-1)*dsmult 736 737 endif 738 goto 300 739 endif 740 enddo 741 742 300 continue 743 744 enddo 745 enddo 746 enddo 747 enddo 748 enddo 749 enddo 750 751 if (.not.ma_pop_stack(l_r3)) 752 1 call errquit('tce_mrcc_iface_r3: MA problem',1,MA_ERR) 753 754 END IF 755 END IF 756 END IF 757 END IF 758 END DO 759 END DO 760 END DO 761 END DO 762 END DO 763 END DO 764 endif 765c 766c --------------- 767c R4 active 768c --------------- 769c 770 if(needr4act) then 771c DO p5b = noab+1,noab+nvab 772c DO p6b = p5b,noab+nvab 773c DO p7b = noab+1,noab+nvab 774c DO p8b = p7b,noab+nvab 775c DO h1b = 1,noab 776c DO h2b = 1,noab 777c DO h3b = h2b,noab 778c DO h4b = h3b,noab 779 780 DO p5b = noab+1,noab+nvab 781 DO p6b = p5b,noab+nvab 782 DO p7b = p6b,noab+nvab 783 DO p8b = p7b,noab+nvab 784 DO h1b = 1,noab 785 DO h2b = h1b,noab 786 DO h3b = h2b,noab 787 DO h4b = h3b,noab 788 789 IF ((.not.restricted).or.(int_mb(k_spin+p5b-1)+int_mb(k_spin+p6b-1 790 &)+int_mb(k_spin+p7b-1)+int_mb(k_spin+p8b-1)+int_mb(k_spin+h1b-1)+i 791 &nt_mb(k_spin+h2b-1)+int_mb(k_spin+h3b-1)+int_mb(k_spin+h4b-1).ne.1 792 &6)) THEN 793 IF (int_mb(k_spin+p5b-1)+int_mb(k_spin+p6b-1)+int_mb(k_spin+p7b-1) 794 &+int_mb(k_spin+p8b-1) .eq. int_mb(k_spin+h1b-1)+int_mb(k_spin+h2b- 795 &1)+int_mb(k_spin+h3b-1)+int_mb(k_spin+h4b-1)) THEN 796 IF (ieor(int_mb(k_sym+p5b-1),ieor(int_mb(k_sym+p6b-1),ieor(int_mb( 797 &k_sym+p7b-1),ieor(int_mb(k_sym+p8b-1),ieor(int_mb(k_sym+h1b-1),ieo 798 &r(int_mb(k_sym+h2b-1),ieor(int_mb(k_sym+h3b-1),int_mb(k_sym+h4b-1) 799 &))))))) .eq. ieor(irrep_v,ieor(irrep_t,irrep_t))) THEN 800 IF ((log_mb(k_active+p5b-1).eqv..true.).and.(log_mb(k_active+p6b-1 801 &).eqv..true.).and.(log_mb(k_active+p7b-1).eqv..true.).and.(log_mb( 802 &k_active+p8b-1).eqv..true.).and.(log_mb(k_active+h1b-1).eqv..true. 803 &).and.(log_mb(k_active+h2b-1).eqv..true.).and.(log_mb(k_active+h3b 804 &-1).eqv..true.).and.(log_mb(k_active+h4b-1).eqv..true.)) THEN 805 806 size = int_mb(k_range+p5b-1) * int_mb(k_range+p6b-1) * int_ 807 &mb(k_range+p7b-1) * int_mb(k_range+p8b-1) * int_mb(k_range+h1b-1) 808 &* int_mb(k_range+h2b-1) * int_mb(k_range+h3b-1) * int_mb(k_range+h 809 &4b-1) 810 811 if (.not.ma_push_get(mt_dbl,size,'r4mi',l_r4,k_r4)) 812 1 call errquit('tce_mrcc_iface_r4: MA problem',0,MA_ERR) 813 814 call get_hash_block(d_r4m(iref),dbl_mb(k_r4),size, 815 1 int_mb(k_r4_offsetm(iref)),(h4b - 1 + noab * (h3b - 1 + noab * 816 1(h2b - 1 + noab * (h1b - 1 + noab * (p8b - noab - 1 + nvab * (p7b 817 1 - noab - 1 + nvab * (p6b - noab - 1 + nvab * (p5b - noab - 1)))) 818 1))))) 819 820 m = 0 821 822 orbspin(1) = int_mb(k_spin+p5b-1)-1 823 orbspin(2) = int_mb(k_spin+p6b-1)-1 824 orbspin(3) = int_mb(k_spin+p7b-1)-1 825 orbspin(4) = int_mb(k_spin+p8b-1)-1 826 orbspin(5) = int_mb(k_spin+h1b-1)-1 827 orbspin(6) = int_mb(k_spin+h2b-1)-1 828 orbspin(7) = int_mb(k_spin+h3b-1)-1 829 orbspin(8) = int_mb(k_spin+h4b-1)-1 830 831 do p5=1,int_mb(k_range+p5b-1) 832 do p6=1,int_mb(k_range+p6b-1) 833 do p7=1,int_mb(k_range+p7b-1) 834 do p8=1,int_mb(k_range+p8b-1) 835 do h1=1,int_mb(k_range+h1b-1) 836 do h2=1,int_mb(k_range+h2b-1) 837 do h3=1,int_mb(k_range+h3b-1) 838 do h4=1,int_mb(k_range+h4b-1) 839 840 m = m + 1 841 842 orbindex(1) = (1 - orbspin(1)+ 843 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p5b-1)+p5-1))/2 844 orbindex(2) = (1 - orbspin(2)+ 845 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p6b-1)+p6-1))/2 846 orbindex(3) = (1 - orbspin(3)+ 847 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p7b-1)+p7-1))/2 848 orbindex(4) = (1 - orbspin(4)+ 849 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p8b-1)+p8-1))/2 850 orbindex(5) = (1 - orbspin(5)+ 851 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h1b-1)+h1-1))/2 852 orbindex(6) = (1 - orbspin(6)+ 853 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h2b-1)+h2-1))/2 854 orbindex(7) = (1 - orbspin(7)+ 855 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h3b-1)+h3-1))/2 856 orbindex(8) = (1 - orbspin(8)+ 857 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h4b-1)+h4-1))/2 858 859 orbindex(1) = moindexes(orbindex(1),orbspin(1)+1,iref) 860 orbindex(2) = moindexes(orbindex(2),orbspin(2)+1,iref) 861 orbindex(3) = moindexes(orbindex(3),orbspin(3)+1,iref) 862 orbindex(4) = moindexes(orbindex(4),orbspin(4)+1,iref) 863 orbindex(5) = moindexes(orbindex(5),orbspin(5)+1,iref) 864 orbindex(6) = moindexes(orbindex(6),orbspin(6)+1,iref) 865 orbindex(7) = moindexes(orbindex(7),orbspin(7)+1,iref) 866 orbindex(8) = moindexes(orbindex(8),orbspin(8)+1,iref) 867 868 do iu=1,2 869 do n=1,nbf 870 ioccnew(n,iu) = iocc(n,iref,iu) 871 enddo 872 enddo 873 874 if((iocc(orbindex(1),iref,orbspin(1)+1).eq. 875 1 iocc(orbindex(5),iref,orbspin(5)+1)).or. 876 2 (iocc(orbindex(2),iref,orbspin(2)+1).eq. 877 3 iocc(orbindex(6),iref,orbspin(6)+1)).or. 878 4 (iocc(orbindex(3),iref,orbspin(3)+1).eq. 879 1 iocc(orbindex(7),iref,orbspin(7)+1)).or. 880 2 (iocc(orbindex(4),iref,orbspin(4)+1).eq. 881 3 iocc(orbindex(8),iref,orbspin(8)+1))) then 882 goto 400 883 endif 884 885 if((orbspin(1).ne.orbspin(5)).or. 886 1 (orbspin(2).ne.orbspin(6)).or. 887 2 (orbspin(3).ne.orbspin(7)).or. 888 3 (orbspin(4).ne.orbspin(8))) then 889 goto 400 890 endif 891 892 if( 893 1 ((orbindex(1).eq.orbindex(2)).and.(orbspin(1).eq.orbspin(2))).or. 894 1 ((orbindex(1).eq.orbindex(3)).and.(orbspin(1).eq.orbspin(3))).or. 895 1 ((orbindex(1).eq.orbindex(4)).and.(orbspin(1).eq.orbspin(4))).or. 896 1 ((orbindex(2).eq.orbindex(3)).and.(orbspin(2).eq.orbspin(3))).or. 897 1 ((orbindex(2).eq.orbindex(4)).and.(orbspin(2).eq.orbspin(4))).or. 898 1 ((orbindex(3).eq.orbindex(4)).and.(orbspin(3).eq.orbspin(4))).or. 899 1 ((orbindex(5).eq.orbindex(6)).and.(orbspin(5).eq.orbspin(6))).or. 900 1 ((orbindex(5).eq.orbindex(7)).and.(orbspin(5).eq.orbspin(7))).or. 901 1 ((orbindex(5).eq.orbindex(8)).and.(orbspin(5).eq.orbspin(8))).or. 902 1 ((orbindex(6).eq.orbindex(7)).and.(orbspin(6).eq.orbspin(7))).or. 903 1 ((orbindex(6).eq.orbindex(8)).and.(orbspin(6).eq.orbspin(8))).or. 904 1 ((orbindex(7).eq.orbindex(8)).and.(orbspin(7).eq.orbspin(8))) 905 1 ) then 906 goto 400 907 endif 908 909 ioccnew(orbindex(1),orbspin(1)+1) = iocc(orbindex(5),iref, 910 1 orbspin(5)+1) 911 ioccnew(orbindex(2),orbspin(2)+1) = iocc(orbindex(6),iref, 912 1 orbspin(6)+1) 913 ioccnew(orbindex(3),orbspin(3)+1) = iocc(orbindex(7),iref, 914 1 orbspin(7)+1) 915 ioccnew(orbindex(4),orbspin(4)+1) = iocc(orbindex(8),iref, 916 1 orbspin(8)+1) 917 ioccnew(orbindex(5),orbspin(5)+1) = iocc(orbindex(1),iref, 918 1 orbspin(1)+1) 919 ioccnew(orbindex(6),orbspin(6)+1) = iocc(orbindex(2),iref, 920 1 orbspin(2)+1) 921 ioccnew(orbindex(7),orbspin(7)+1) = iocc(orbindex(3),iref, 922 1 orbspin(3)+1) 923 ioccnew(orbindex(8),orbspin(8)+1) = iocc(orbindex(4),iref, 924 1 orbspin(4)+1) 925 926 do n=1,nref 927 isfound = .true. 928 do iu=1,2 929 do o=1,nbf 930 isfound = isfound.and.(iocc(o,n,iu).eq.ioccnew(o,iu)) 931 enddo 932 enddo 933 if(isfound) then 934c write(LuOut,"('Internal amplitude',I4,'->',I4,2F16.12)")iref,n, 935c 1 dbl_mb(k_r4+m-1) 936 if(iref.ne.n) then 937c =============== 938 dist = 0 939 do iu=1,2 940 do i1=1,nmo(iu) 941 ioccnew(i1,iu) = iocc(i1,iref,iu) 942 enddo 943 enddo 944 945 i2 = 0 946 do i1=min(orbindex(1),orbindex(5)), 947 1 max(orbindex(1),orbindex(5)) 948 i2 = i2 + 1 949 if(i2.lt.abs(orbindex(1)-orbindex(5))) then 950 if(iocc(i1+1,iref,orbspin(1)+1).eq.1) then 951 dist=dist+1 952 endif 953 endif 954 enddo 955 956 ioccnew(orbindex(1),orbspin(1)+1) = iocc(orbindex(5),iref, 957 1 orbspin(5)+1) 958 ioccnew(orbindex(5),orbspin(5)+1) = iocc(orbindex(1),iref, 959 1 orbspin(1)+1) 960 961 i2 = 0 962 do i1=min(orbindex(2),orbindex(6)), 963 1 max(orbindex(2),orbindex(6)) 964 i2 = i2 + 1 965 if(i2.lt.abs(orbindex(2)-orbindex(6))) then 966 if(iocc(i1+1,iref,orbspin(2)+1).eq.1) then 967 dist=dist+1 968 endif 969 endif 970 enddo 971 972 ioccnew(orbindex(2),orbspin(2)+1) = iocc(orbindex(6),iref, 973 1 orbspin(6)+1) 974 ioccnew(orbindex(6),orbspin(6)+1) = iocc(orbindex(2),iref, 975 1 orbspin(2)+1) 976 977 i2 = 0 978 do i1=min(orbindex(3),orbindex(7)), 979 1 max(orbindex(3),orbindex(7)) 980 i2 = i2 + 1 981 if(i2.lt.abs(orbindex(3)-orbindex(7))) then 982 if(iocc(i1+1,iref,orbspin(3)+1).eq.1) then 983 dist=dist+1 984 endif 985 endif 986 enddo 987 988 ioccnew(orbindex(3),orbspin(3)+1) = iocc(orbindex(7),iref, 989 1 orbspin(7)+1) 990 ioccnew(orbindex(7),orbspin(7)+1) = iocc(orbindex(3),iref, 991 1 orbspin(3)+1) 992 993 i2 = 0 994 do i1=min(orbindex(4),orbindex(8)), 995 1 max(orbindex(4),orbindex(8)) 996 i2 = i2 + 1 997 if(i2.lt.abs(orbindex(4)-orbindex(8))) then 998 if(iocc(i1+1,iref,orbspin(4)+1).eq.1) then 999 dist=dist+1 1000 endif 1001 endif 1002 enddo 1003 1004 dsmult = 1.0d0 1005 if(mod(dist,2).eq.1) dsmult = -dsmult 1006 1007c if(nodezero)then 1008c write(6,"('T4 iref/n:',2I4,f4.1)")iref,n,dsmult 1009c endif 1010c =============== 1011 dbl_mb(k_heff+n-1+(iref-1)*nref) = 1012 1 dbl_mb(k_r4+m-1)*dsmult 1013 endif 1014 goto 400 1015 endif 1016 enddo 1017 1018 400 continue 1019 1020 enddo 1021 enddo 1022 enddo 1023 enddo 1024 enddo 1025 enddo 1026 enddo 1027 enddo 1028 1029 if (.not.ma_pop_stack(l_r4)) 1030 1 call errquit('tce_mrcc_iface_r4: MA problem',1,MA_ERR) 1031 1032 END IF 1033 END IF 1034 END IF 1035 END IF 1036 END DO 1037 END DO 1038 END DO 1039 END DO 1040 END DO 1041 END DO 1042 END DO 1043 END DO 1044 1045 endif 1046 1047 endif !sub 1048 1049 enddo !iref 1050 1051 return 1052 end 1053c 1054c ============================================== 1055c Diagonalize effective Hamiltonian 1056c ============================================== 1057c 1058 subroutine tce_diagonalize_heff(rtdb,iter) 1059 implicit none 1060#include "global.fh" 1061#include "mafdecls.fh" 1062#include "sym.fh" 1063#include "util.fh" 1064#include "stdio.fh" 1065#include "errquit.fh" 1066#include "tce.fh" 1067#include "tce_mrcc.fh" 1068#include "tce_main.fh" 1069#include "rtdb.fh" 1070#include "tcgmsg.fh" 1071#include "msgtypesf.h" 1072#include "msgids.fh" 1073 1074c integer nref 1075 double precision heff(nref,nref) 1076 double precision vl(nref,nref) 1077 double precision vr(nref,nref) 1078 double precision er(nref) 1079 double precision ei(nref) 1080 double precision work(4*nref) 1081 integer info 1082 integer i,j,l 1083 double precision x 1084 integer k 1085 logical nodezero 1086 integer rtdb,itarget 1087c integer nrootmuc 1088 double precision dvalue,dsum 1089 double precision dfin 1090 character*4 ds 1091 integer l_buff,k_buff 1092c integer iignore 1093 double precision isum 1094 integer ddblsize,inntsize 1095 1096ckbn rootfromoverlap 1097 double precision rootoverlap,rootoverlapmax 1098 integer stateselectedbyoverlap 1099 integer iter 1100 1101 nodezero = (ga_nodeid().eq.0) 1102 ddblsize=MA_sizeof(MT_DBL,1,MT_BYTE) 1103 inntsize=MA_sizeof(MT_INT,1,MT_BYTE) 1104 1105 call ga_sync() 1106 1107 if(lusesub) then 1108 do i=1,nref*nref 1109 dbl_mb(k_heff+i-1) = 0.0d0 1110 enddo 1111 call ga_sync() 1112 call ga_get(g_heff,1,nref*nref,1,1, 1113 1 dbl_mb(k_heff),1) 1114c call ga_print(g_heff) 1115 endif 1116 1117 do i=1,nref 1118 do j=1,nref 1119 x = dbl_mb(k_heff+j-1+(i-1)*nref) 1120 heff(j,i) = x 1121 end do 1122 end do 1123 1124 do i=1,nref 1125 er(i) = 0.0d0 1126 ei(i) = 0.0d0 1127 do j=1,nref 1128 vl(i,j)=0.0d0 1129 vr(i,j)=0.0d0 1130 enddo 1131 enddo 1132 1133 if(nodezero.and.(nref.lt.21)) then 1134c call ma_print(dbl_mb(k_heff),nref,nref,'Heff') 1135 1136 write(LuOut,"(/,'Heff',/, 1137 1 '=============================================')") 1138 do i=1,nref 1139 write(LuOut,"(i5,i5,100F14.8)")ga_nodeid(),i, 1140 1 (dbl_mb(k_heff+(j-1)*nref+i-1),j=1,nref) 1141 enddo 1142 1143 endif 1144c call ga_sync() 1145c call util_flush(LuOut) 1146c if((ga_nodeid().eq.5).and.(nref.lt.21)) then 1147c call ma_print(dbl_mb(k_heff),nref,nref,'Heff') 1148c 1149c write(LuOut,"(/,'Heff 5',/, 1150c 1 '=============================================')") 1151c do i=1,nref 1152c write(LuOut,"(i5,i5,100F14.8)")ga_nodeid(),i, 1153c 1 (dbl_mb(k_heff+(j-1)*nref+i-1),j=1,nref) 1154c enddo 1155 1156c endif 1157c call util_flush(LuOut) 1158c call ga_sync() 1159 1160c call util_flush(LuOut) 1161c write(6,*)'BEFORE',ga_nodeid() 1162c call util_flush(LuOut) 1163 call ga_sync() 1164c if(nodezero)write(6,*)'TEST 3' 1165 if(nodezero) then 1166c call DGEEV('V','V',nref,heff,nref,er,ei,vl,nref,vr, 1167 call util_dgeev('V','V',nref,heff,nref,er,ei,vl,nref,vr, 1168 $ nref,work,4*nref,info) 1169 if(info .ne. 0) call errquit('Heff diagonalization',0,CALC_ERR) 1170 call amp_stabilization(vl,vr,nref) 1171 endif 1172c if(nodezero)write(6,*)'TEST 4' 1173 call ga_sync() 1174 1175c call util_flush(LuOut) 1176c write(6,*)'AFTER',ga_nodeid() 1177c call util_flush(LuOut) 1178 1179c if(nodezero.and..not.lconverged) then 1180c call ma_print(dbl_mb(k_heff),nref,nref,'Heff') 1181c endif 1182 1183c if(lconverged.and.nodezero) then 1184 if(nodezero) then 1185 if(nref.lt.21) then 1186 1187 write(LuOut,"(/,'Eigenvalues (real and imaginary)',/, 1188 1 '=============================================')") 1189 do i=1,nref 1190 write(LuOut,"(F18.12,100F14.8)")er(i),ei(i) 1191 enddo 1192 1193 write(LuOut,"(/,'Left eigenvectors',/, 1194 1 '=============================================')") 1195 do i=1,nref 1196 write(LuOut,"(i5,100F14.8)")i,(vl(i,j),j=1,nref) 1197 enddo 1198 1199c write(LuOut,"(/,'Left eigenvectors - squares',/, 1200c 1 '=============================================')") 1201c do i=1,nref 1202c write(LuOut,"(i5,35f18.12)")i,((vl(i,j)*vl(i,j)),j=1,nref) 1203c enddo 1204 1205c endif 1206c endif 1207 1208 write(LuOut,"(/,'Right eigenvectors',/, 1209 1 '=============================================')") 1210 do i=1,nref 1211 write(LuOut,"('VR',i5,100F14.8)")i,(vr(i,j),j=1,nref) 1212 enddo 1213 1214 endif 1215 endif 1216 1217 do i=1,nref 1218 do j=1,nref 1219 dbl_mb(k_sqc+(i-1)*nref+j-1)=vr(i,j) 1220 dbl_mb(k_sqcl+(i-1)*nref+j-1)=vl(i,j) 1221 enddo 1222 enddo 1223 1224 call ga_brdcst(Msg_Vec_EVal+21,dbl_mb(k_sqc), 1225 1 ddblsize*nref*nref, 0) 1226 call ga_sync() 1227 call ga_brdcst(Msg_Vec_EVal+20,dbl_mb(k_sqcl), 1228 1 ddblsize*nref*nref, 0) 1229 call ga_sync() 1230 1231 if (.not.rtdb_get(rtdb,'mrcc:zignore',mt_int,1,iignore)) 1232 1 iignore = 0 1233 1234 if(.not.nodezero) then 1235 do i=1,nref 1236 do j=1,nref 1237 vr(i,j) = dbl_mb(k_sqc+(i-1)*nref+j-1) 1238 vl(i,j) = dbl_mb(k_sqcl+(i-1)*nref+j-1) 1239 enddo 1240 enddo 1241 endif 1242 1243c-------- 1244 1245 epsilon = 0.0d0 1246 do i=1,nref 1247 isum = 0.0d0 1248 do j=1,nref 1249 isum = isum + vr(j,i) 1250 enddo 1251 if((iignore.eq.0).or.(abs(isum).gt.1.0d-5)) then 1252 if(epsilon.gt.min(epsilon,er(i)))mkroot=i 1253 epsilon = min(epsilon,er(i)) 1254 endif 1255 enddo 1256 1257c-------- 1258 1259 if (.not.rtdb_get(rtdb,'mrcc:rootmuc',mt_int,1,nrootmuc)) 1260 1 nrootmuc = 0 1261c ------ 1262 if(nrootmuc.gt.0) then ! 1 1263 dfin = 0.0d0 1264 1265 do j=1,nref 1266 dsum = 0.0d0 1267 isum = 0.0d0 1268 do i=1,nref 1269 isum = isum + vr(i,j) 1270 enddo 1271 if((iignore.eq.0).or.(abs(isum).gt.1.0d-5)) then 1272 do i=1,nref 1273 write(ds,"(I3.3)")i 1274 if (.not.rtdb_get(rtdb,'mrcc:rootmuc'//ds,mt_dbl,1,dvalue)) 1275 1 dvalue = 0.0d0 1276 if(dvalue.lt.0.0d0) then 1277 if(abs(vr(i,j)).lt.abs(dvalue)) 1278 1 dsum = dsum + (dvalue)*(vr(i,j)) 1279 else 1280 dsum = dsum + (dvalue)*(vr(i,j)) 1281 endif 1282 enddo 1283c write(6,"('SUM:',I4,4F16.12)")j, 1284c 1 abs(dsum) 1285 if(dfin.lt.abs(dsum)) then 1286c write(6,"('I am watching reference #',I4,4F16.12)")j 1287 mkroot=j 1288 epsilon = er(j) 1289 dfin = abs(dsum) 1290 endif 1291 endif 1292 enddo 1293 1294 else ! 1 1295 1296 if (rtdb_get(rtdb,'bwcc:targetroot',mt_int,1,itarget)) then 1297c mkroot = itarget 1298 do i=1,nref 1299 1300 isum = 0.0d0 1301 do j=1,nref 1302c if(abs(vr(j,i)).lt.1.0d-8) 1303 isum = isum + vr(j,i) 1304 enddo 1305 if((iignore.eq.0).or.(abs(isum).gt.1.0d-5)) then 1306 1307 dvalue = er(i) 1308 k = 0 1309 do j=1,nref 1310 if(j.ne.i) then 1311 1312 isum = 0.0d0 1313 do l=1,nref 1314c if(abs(vr(l,j)).lt.1.0d-8) 1315 isum = isum + vr(l,j) 1316 enddo 1317 if((iignore.eq.0).or.(abs(isum).gt.1.0d-5)) then 1318 1319 1320 if(er(j).lt.er(i))k=k+1 1321 1322 endif 1323 1324 endif 1325 enddo 1326 if(k.eq.(itarget-1)) then 1327 itarget = i 1328 goto 63454 1329 endif 1330 1331 endif 1332 1333 enddo 1334 133563454 continue 1336 1337 mkroot = itarget 1338 epsilon = er(itarget) 1339 1340 endif 1341 endif 1342 1343 1344ckbn find the requested root according to overlap instead 1345 call ga_sync() 1346 if(iter .gt. iroottooverlapiter) then 1347c if(.not.lconverged) then 1348 if(lrootfromoverlap) then 1349 stateselectedbyoverlap=0 1350 rootoverlapmax=-100.0d0 1351 do i=1,nref 1352 rootoverlap=0.0d0 1353 isum = 0.0d0 1354 do j=1,nref 1355c if(nodezero)write(*,'(A,2I5,A,F17.10)')"vr(",j,i,")",vr(j,i) 1356 rootoverlap=rootoverlap+bwcoefwanted(j)*vr(j,i) 1357 isum = isum + vr(j,i) 1358 enddo 1359 rootoverlap=abs(rootoverlap) 1360c write(*,'(A,F17.10)') "Overlap ", rootovelap 1361c write(*,*) "i1 ",i," rootoverlapmax ",rootoverlapmax,stateselectedbyoverlap,i 1362 if(nodezero) write(*,'(A,I5,F17.10,A,F17.10)') "Overlap ",i, 1363 + rootoverlap, " isum ",isum 1364 if((iignore.eq.0).or.(abs(isum).gt.1.0d-5)) then 1365 if(rootoverlap .gt. rootoverlapmax) then 1366 rootoverlapmax=rootoverlap 1367 stateselectedbyoverlap=i 1368c write(*,*) "i2 ",i," rootovelapmax ",rootovelapmax,stateselectedbyoverlap,i 1369 endif 1370 endif 1371 enddo 1372 mkroot= stateselectedbyoverlap 1373 epsilon = er(stateselectedbyoverlap) 1374 if(nodezero) write(LuOut,'(A,I5)') 1375 + 'MRCC root selected by overlap is ',stateselectedbyoverlap 1376 endif ! lrootfromoverlap 1377c else ! lconverged 1378c if(nodezero) write(*,'(A,L3,A,I5)') "Converged ",lconverged, 1379c + " so keeping old root ", mkrootold 1380c endif !lconverged 1381 endif 1382 call ga_sync() 1383 1384 1385 1386 1387 1388 1389 if (nodezero.and.(nref.lt.21)) then 1390 write(6,"('Target root: ',I4)")mkroot 1391 end if 1392 1393ckbn Introduce some checks before proceeding further. 1394 1395ckbn check of mkroot, crash if mkroot gt nref 1396 if(nodezero) then 1397 if(mkroot.gt.nref) call errquit 1398 + ('root followed is greater than total number of references',0, 1399 + CALC_ERR) 1400 endif 1401 1402 1403ckbn check whether there is imaginary eigen values 1404 do i=1,nref 1405c write(LuOut,'(A,F17.10)')'Imaginary eigenvalue',ei(i) 1406c if (nodezero) call util_flush(LuOut) 1407 if(abs(ei(i)).ge.toleiimag) then 1408 write(LuOut,'(A,F15.10)') 1409 + 'Warning: complex Heff eigenvalue detected',ei(i) 1410 if(i.eq.mkroot) then 1411 write(LuOut,*) "ignorecomplex1 ", ignorecomplex 1412c if (rtdb_get(rtdb,'mrcc:ignorecomplex',mt_log,1, 1413c + ignorecomplex)) ignorecomplex = .true. 1414c call errquit('Complex root found and no ignorecomplex option', 1415c + 0,RTDB_ERR) 1416c else 1417c if(nodezero) write(*,*) "ignorecomplex1. ", ignorecomplex 1418c ignorecomplex = .true. 1419c endif 1420c if(nodezero) write(*,*) "ignorecomplex2 ", ignorecomplex 1421 if(ignorecomplex) then 1422 write(LuOut,'(A,F15.10)') 1423 + 'Warning: Proceeding with complex Heff eigenvalue ',ei(i) 1424 else 1425 call errquit('Warning:complex Heff eigenvalue detected',0, 1426 + UNKNOWN_ERR) 1427 endif 1428 endif 1429 endif 1430 enddo 1431 1432 1433 1434 1435 if(nref.gt.20) then 1436 if(nodezero)write(LuOut,"(/,'Right eigenvector for target',I7, 1437 1 /)")mkroot 1438c if (.not.ma_push_get(mt_dbl,nref,'buff',l_buff,k_buff)) 1439c 1 call errquit('tce_mrcc_iface_buff: MA problem',0,MA_ERR) 1440 do i=1,nref 1441 if(abs(vr(i,mkroot)).gt.0.05d0) then 1442 if(nodezero)write(LuOut,"(I7,' ',F11.8)")i,vr(i,mkroot) 1443 endif 1444 enddo 1445c if (.not.ma_pop_stack(l_buff)) 1446c 1 call errquit('tce_mrcc_iface_buff: MA problem',1,MA_ERR) 1447 1448 endif 1449 1450 if (nodezero) call util_flush(LuOut) 1451 1452cjb broadcasts 1453 1454 call ga_brdcst(Msg_Vec_EVal+MSGINT+30,mkroot,inntsize, 0) 1455 call ga_brdcst(Msg_Vec_EVal+MSGDBL+31,epsilon,ddblsize, 0) 1456 call ga_sync() 1457 1458c write(LuOut,"(/,'Epsilon: ',2F16.12,/)") epsilon 1459 1460 return 1461 end 1462 1463c 1464c ============================================== 1465c Clean internal amplitudes 1466c ============================================== 1467c 1468 subroutine tce_internal_t_zero(d_t1m,d_t2m,k_t1_offsetm, 1469 1 k_t2_offsetm,lneedt3,d_t3m,k_t3_offsetm,rtdb) 1470c 1 k_t2_offsetm,nref,lneedt3,d_t3m,k_t3_offsetm,rtdb) 1471 implicit none 1472#include "global.fh" 1473#include "mafdecls.fh" 1474#include "sym.fh" 1475#include "util.fh" 1476#include "stdio.fh" 1477#include "errquit.fh" 1478#include "tce.fh" 1479#include "tce_mrcc.fh" 1480#include "tce_main.fh" 1481#include "rtdb.fh" 1482 1483 integer d_t1m(maxref),d_t2m(maxref) 1484 integer d_t3m(maxref) 1485 integer k_t1_offsetm(maxref),k_t2_offsetm(maxref) 1486 integer k_t3_offsetm(maxref) 1487 integer size,p5b,h6b 1488c integer nref 1489 integer l_t1,k_t1 1490 integer l_t2,k_t2 1491 integer l_t3,k_t3 1492 integer i,j,k,l,m 1493 integer iref,rtdb 1494 integer p1b,p2b,h3b,h4b 1495 integer p3b,p4b,h5b,p6b 1496 integer h1b,h2b 1497c logical lneedt3,limprovet 1498 logical lneedt3 1499 integer orbindex(6),orbspin(6) 1500 EXTERNAL NXTASKsub 1501 EXTERNAL NXTASK 1502 INTEGER NXTASKsub 1503 INTEGER NXTASK 1504 INTEGER nxt 1505 INTEGER nprocs 1506 INTEGER count 1507 1508c if (.not.rtdb_get(rtdb,'mrcc:improvetiling',mt_log,1,limprovet)) 1509c 1 limprovet = .false. 1510 1511 do iref=1,nref 1512 1513 if((int_mb(k_refafi+iref-1).eq.int_mb(k_innodes 1514 1 +ga_nnodes()+ga_nodeid())).or.(.not.lusesub)) then 1515 1516 k_sym = k_symm(iref) 1517 k_offset = k_offsetm(iref) 1518 k_range = k_rangem(iref) 1519 k_spin = k_spinm(iref) 1520 k_movecs_sorted = k_movecs_sortedm(iref) 1521 1522 noa = nblcks(1,iref) 1523 nob = nblcks(2,iref) 1524 nva = nblcks(3,iref) 1525 nvb = nblcks(4,iref) 1526 1527 noab = noa+nob 1528 nvab = nva+nvb 1529 1530 do p5b = noab+1,noab+nvab 1531 do h6b = 1,noab 1532 1533 if (int_mb(k_spin+p5b-1) .eq. int_mb(k_spin+h6b-1)) then 1534 if (ieor(int_mb(k_sym+p5b-1),int_mb(k_sym+h6b-1)).eq.irrep_t)then 1535 if ((.not.restricted).or.(int_mb(k_spin+p5b-1)+int_mb(k_spin+h6b-1 1536 &).ne.4)) then 1537 if(log_mb(k_isactive(iref)+p5b-1).and. 1538 &log_mb(k_isactive(iref)+h6b-1)) then 1539 1540 size = int_mb(k_range+h6b-1) * int_mb(k_range+p5b-1) 1541 if (.not.ma_push_get(mt_dbl,size,'t1',l_t1,k_t1)) 1542 1 call errquit('tce_mrcc_iface_t1: MA problem',0,MA_ERR) 1543 1544c call dfill(size,0.0d0,dbl_mb(k_t1),1) 1545 do i=1,size 1546 dbl_mb(k_t1+i-1)=0.0d0 1547 enddo 1548 1549cjb ============================ 1550 1551 if(limprovet) then 1552 1553 call get_hash_block(d_t1m(iref),dbl_mb(k_t1),size, 1554 1 int_mb(k_t1_offsetm(iref)),h6b-1+noab*(p5b-noab-1)) 1555 1556 k=0 1557 do i=1,int_mb(k_range+p5b-1) 1558 do j=1,int_mb(k_range+h6b-1) 1559 k = k + 1 1560 1561 orbspin(1) = int_mb(k_spin+p5b-1)-1 1562 orbspin(2) = int_mb(k_spin+h6b-1)-1 1563 1564 orbindex(1) = (1 - orbspin(1)+ 1565 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p5b-1)+i-1))/2 1566 orbindex(2) = (1 - orbspin(2)+ 1567 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h6b-1)+j-1))/2 1568 1569 orbindex(1) = moindexes(orbindex(1),orbspin(1)+1,iref) 1570 orbindex(2) = moindexes(orbindex(2),orbspin(2)+1,iref) 1571 1572 if(isactive(orbindex(1),orbspin(1)+1).and. 1573 1 isactive(orbindex(2),orbspin(2)+1).or.(.not.limprovet)) then 1574 1575 dbl_mb(k_t1+k-1) = 0.0d0 1576 1577 endif 1578 enddo 1579 enddo 1580 1581 endif ! limprovet 1582 1583c ============================= 1584 1585c do i=1,size 1586c dbl_mb(k_t1+i-1) = 0.0d0 1587c enddo 1588 1589 call put_hash_block(d_t1m(iref),dbl_mb(k_t1),size, 1590 1 int_mb(k_t1_offsetm(iref)),((p5b-noab-1)*noab+h6b-1)) 1591 1592 if (.not.ma_pop_stack(l_t1)) 1593 1 call errquit('tce_mrcc_iface_t1: MA problem',1,MA_ERR) 1594 1595 endif 1596 endif 1597 endif 1598 endif 1599 enddo 1600 enddo 1601 1602cjb new in parallel 1603 1604 nxt = 0 1605 if(limprovet) then 1606 if(lusesub) then 1607 nprocs=GA_pgroup_NNODES(mypgid) 1608 nxt=NXTASKsub(nprocs,1,mypgid) 1609 else 1610 nprocs=GA_NNODES() 1611 nxt=NXTASK(nprocs,1) 1612 endif 1613 count = 0 1614 endif 1615 1616 DO p1b = noab+1,noab+nvab 1617 DO p2b = p1b,noab+nvab 1618 DO h3b = 1,noab 1619 DO h4b = h3b,noab 1620 1621 IF ((nxt.eq.count).or.(.not.limprovet)) THEN 1622 1623 IF (int_mb(k_spin+p1b-1)+int_mb(k_spin+p2b-1) .eq. int_mb(k_spin+h 1624 &3b-1)+int_mb(k_spin+h4b-1)) THEN 1625 IF (ieor(int_mb(k_sym+p1b-1),ieor(int_mb(k_sym+p2b-1),ieor(int_mb( 1626 &k_sym+h3b-1),int_mb(k_sym+h4b-1)))) .eq. irrep_t) THEN 1627 IF ((.not.restricted).or.(int_mb(k_spin+p1b-1)+int_mb(k_spin+p2b-1 1628 &)+int_mb(k_spin+h3b-1)+int_mb(k_spin+h4b-1).ne.8)) THEN 1629 if(log_mb(k_isactive(iref)+p1b-1).and. 1630 1 log_mb(k_isactive(iref)+p2b-1).and. 1631 2 log_mb(k_isactive(iref)+h3b-1).and. 1632 3 log_mb(k_isactive(iref)+h4b-1)) then 1633 1634 size = int_mb(k_range+p1b-1) * int_mb(k_range+p2b-1) * int_ 1635 &mb(k_range+h3b-1) * int_mb(k_range+h4b-1) 1636 1637 if (.not.ma_push_get(mt_dbl,size,'t2',l_t2,k_t2)) 1638 1 call errquit('tce_mrcc_iface_t2: MA problem',0,MA_ERR) 1639 1640c call dfill(size,0.0d0,dbl_mb(k_t2),1) 1641 do i=1,size 1642 dbl_mb(k_t2+i-1)=0.0d0 1643 enddo 1644 1645c =============================================================== 1646 1647 if(limprovet) then 1648 1649 call get_hash_block(d_t2m(iref),dbl_mb(k_t2),size, 1650 1 int_mb(k_t2_offsetm(iref)),h4b-1+noab*(h3b-1+noab*(p2b- 1651 &noab-1+nvab*(p1b - noab - 1)))) 1652 1653 m = 0 1654 do i=1,int_mb(k_range+p1b-1) 1655 do j=1,int_mb(k_range+p2b-1) 1656 do k=1,int_mb(k_range+h3b-1) 1657 do l=1,int_mb(k_range+h4b-1) 1658 m = m + 1 1659 1660 orbspin(1) = int_mb(k_spin+p1b-1)-1 1661 orbspin(2) = int_mb(k_spin+p2b-1)-1 1662 orbspin(3) = int_mb(k_spin+h3b-1)-1 1663 orbspin(4) = int_mb(k_spin+h4b-1)-1 1664 1665 orbindex(1) = (1 - orbspin(1)+ 1666 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p1b-1)+i-1))/2 1667 orbindex(2) = (1 - orbspin(2)+ 1668 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+p2b-1)+j-1))/2 1669 orbindex(3) = (1 - orbspin(3)+ 1670 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h3b-1)+k-1))/2 1671 orbindex(4) = (1 - orbspin(4)+ 1672 1 int_mb(k_mo_indexm(iref)+int_mb(k_offset+h4b-1)+l-1))/2 1673 1674 orbindex(1) = moindexes(orbindex(1),orbspin(1)+1,iref) 1675 orbindex(2) = moindexes(orbindex(2),orbspin(2)+1,iref) 1676 orbindex(3) = moindexes(orbindex(3),orbspin(3)+1,iref) 1677 orbindex(4) = moindexes(orbindex(4),orbspin(4)+1,iref) 1678 1679 if(isactive(orbindex(1),orbspin(1)+1).and. 1680 1 isactive(orbindex(2),orbspin(2)+1).and. 1681 2 isactive(orbindex(3),orbspin(3)+1).and. 1682 3 isactive(orbindex(4),orbspin(4)+1).or.(.not.limprovet)) then 1683 1684 dbl_mb(k_t2+m-1)=0.0d0 1685 1686 endif 1687 enddo 1688 enddo 1689 enddo 1690 enddo 1691 1692 endif 1693 1694c =============================================================== 1695 1696c do i=1,size 1697c write(LuOut,*)dbl_mb(k_t2+i-1),'->','0.00000000' 1698c dbl_mb(k_t2+i-1) = 0.0d0 1699c enddo 1700 1701 call put_hash_block(d_t2m(iref),dbl_mb(k_t2),size, 1702 1 int_mb(k_t2_offsetm(iref)),((((p1b-noab-1)*nvab+p2b-noab-1) 1703 2 *noab+h3b-1)*noab+h4b-1)) 1704 1705 if (.not.ma_pop_stack(l_t2)) 1706 1 call errquit('tce_mrcc_iface_t2: MA problem',1,MA_ERR) 1707 1708 END IF 1709 END IF 1710 END IF 1711 endif 1712 1713 if(limprovet) then 1714 if(lusesub) then 1715 nxt=NXTASKsub(nprocs,1,mypgid) 1716 else 1717 nxt=NXTASK(nprocs,1) 1718 endif 1719 endif 1720 1721 endif 1722 1723 if(limprovet)count = count + 1 1724 1725 END DO 1726 END DO 1727 END DO 1728 END DO 1729 1730 if(limprovet) then 1731 if(lusesub) then 1732 nxt=NXTASKsub(-nprocs,1,mypgid) 1733 call GA_Pgroup_SYNC(mypgid) 1734 else 1735 nxt=NXTASK(-nprocs,1) 1736 call GA_SYNC() 1737 endif 1738 endif 1739 1740 if(lneedt3) then 1741 1742 DO p2b = noab+1,noab+nvab 1743 DO p3b = p2b,noab+nvab 1744 DO p4b = p3b,noab+nvab 1745 DO h1b = 1,noab 1746 DO h5b = h1b,noab 1747 DO h6b = h5b,noab 1748 IF (int_mb(k_spin+p2b-1)+int_mb(k_spin+p3b-1)+int_mb(k_spin+p4b-1) 1749 & .eq. int_mb(k_spin+h1b-1)+int_mb(k_spin+h5b-1)+int_mb(k_spin+h6b- 1750 &1)) THEN 1751 IF (ieor(int_mb(k_sym+p2b-1),ieor(int_mb(k_sym+p3b-1),ieor(int_mb( 1752 &k_sym+p4b-1),ieor(int_mb(k_sym+h1b-1),ieor(int_mb(k_sym+h5b-1),int 1753 &_mb(k_sym+h6b-1)))))) .eq. irrep_t) THEN 1754 IF ((.not.restricted).or.(int_mb(k_spin+p2b-1)+int_mb(k_spin+p3b-1 1755 &)+int_mb(k_spin+p4b-1)+int_mb(k_spin+h1b-1)+int_mb(k_spin+h5b-1)+i 1756 &nt_mb(k_spin+h6b-1).ne.12)) THEN 1757 IF ((log_mb(k_isactive(iref)+p2b-1).eqv..true.).and. 1758 1 (log_mb(k_isactive(iref)+p3b-1).eqv..true.).and. 1759 2 (log_mb(k_isactive(iref)+p4b-1).eqv..true.).and. 1760 3 (log_mb(k_isactive(iref)+h1b-1).eqv..true.).and. 1761 4 (log_mb(k_isactive(iref)+h5b-1).eqv..true.).and. 1762 5 (log_mb(k_isactive(iref)+h6b-1).eqv..true.)) THEN 1763 1764 size = int_mb(k_range+p2b-1) * int_mb(k_range+p3b-1) * int_ 1765 &mb(k_range+p4b-1) * int_mb(k_range+h1b-1) * int_mb(k_range+h5b-1) 1766 &* int_mb(k_range+h6b-1) 1767 1768 if (.not.ma_push_get(mt_dbl,size,'t3',l_t3,k_t3)) 1769 1 call errquit('tce_mrcc_iface_t3: MA problem',0,MA_ERR) 1770 1771 do i=1,size 1772 dbl_mb(k_t3+i-1) = 0.0d0 1773 enddo 1774 call put_hash_block(d_t3m(iref),dbl_mb(k_t3),size, 1775 1 int_mb(k_t3_offsetm(iref)),(h6b - 1 + noab * 1776 2 (h5b - 1 + noab * (h1b 1777 &- 1 + noab * (p4b - noab - 1 + nvab * (p3b - noab - 1 + nvab * (p2 1778 &b - noab - 1))))))) 1779 1780 if (.not.ma_pop_stack(l_t3)) 1781 1 call errquit('tce_mrcc_iface_t3: MA problem',1,MA_ERR) 1782 1783 END IF 1784 END IF 1785 END IF 1786 END IF 1787 END DO 1788 END DO 1789 END DO 1790 END DO 1791 END DO 1792 END DO 1793 1794 endif !needt3 1795 1796c write(6,"('CPU BEFORE',I4,F16.12)")ga_nodeid(),util_cpusec() 1797c call ga_sync() 1798 if(lusesub) call ga_pgroup_sync( 1799 1 int_mb(k_innodes+ga_nnodes()+ga_nodeid())) 1800c write(6,"('CPU AFTER',I4,F16.12)")ga_nodeid(),util_cpusec() 1801 1802 endif !sub 1803 1804 enddo !iref 1805 1806 return 1807 end 1808 1809 1810c $Id$ 1811