1 SUBROUTINE tce_ao2e_cholesky_orb(rtdb,d_v2,kax_v2_alpha_offset, 2 1 size_2e) 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 t ( p1 p2 h3 h4 )_t 7 IMPLICIT NONE 8#include "rtdb.fh" 9#include "global.fh" 10#include "mafdecls.fh" 11#include "stdio.fh" 12#include "util.fh" 13#include "bas.fh" 14#include "schwarz.fh" 15#include "sym.fh" 16#include "sf.fh" 17#include "errquit.fh" 18#include "tce.fh" 19#include "tce_main.fh" 20c 21c 22 integer rtdb ! Run-time database 23 integer d_v2 ! MO integrals 24 integer kax_v2_alpha_offset ! MO integrals offset 25 integer size_2e ! 2e file size 26c 27 INTEGER size_2g2a,l_2g2a,k_2g2a 28 INTEGER azone1,azone2,azone3,azone4 29 INTEGER g1b,g2b,g3b,g4b 30 INTEGER igi1,igi2,igi3,igi4 31 INTEGER ii,i,j,k,l,N,ipos1,ipos2,ipos3,ipos4 32 INTEGER del1,del2,p1rel,p2rel 33 INTEGER size_4a,l_4a,k_4a 34c 35 integer mu,nu,rho,sigma 36 integer mu_lo,mu_hi 37 integer nu_lo,nu_hi 38 integer rho_lo,rho_hi 39 integer sigma_lo,sigma_hi 40 integer mu_range 41 integer nu_range 42 integer rho_range 43 integer sigma_range 44 integer mu1,nu1,rho1,sigma1 45 integer shift_mu,shift_nu 46 integer shift_rho,shift_sigma 47 integer work1,work2 ! Work array sizes 48 integer l_work1,k_work1 ! Work array 1 49 integer l_work2,k_work2 ! Work array 2 50 integer imu1,inu1,irho1,isigma1 51c 52 integer l_movecs_orb,k_movecs_orb 53 integer l_gpair,k_gpair 54 integer len_pair,g12_shift 55c ATTENTION,ACHTUNG,UWAGA 2000 - max # of CPU 56c 57 integer size_2g2z,l_2g2z,k_2g2z 58 integer tot_azone1_sh,tot_azone2_sh 59 integer tot_azone3_sh,tot_azone4_sh 60 integer ixi,jxi,point_pair 61c 62 integer iha,ihb !number of corr. alpha, beta holes 63 integer ipa,ipb !number of corr. alpha, beta particles 64c compression 65 integer max_size_temp,size_temp,sumx 66 double precision wall,cpu,wall1,cpu1,wall2,cpu2,wall3,cpu3 67c 68 double precision tot_zone(1000) 69c 70 integer l_3a1m_offset,k_3a1m_offset,size_3a1mf 71 integer l_1a3m_offset,k_1a3m_offset,size_1a3mf 72 integer d_1a3m,d_3a1m 73 integer size_3a1m 74 integer l_3a1m,k_3a1m 75 integer size_amc 76 integer key_3a1m,offset_3a1m 77 integer l_4a_sort,k_4a_sort 78 integer size_2a2m,l_2a2m,k_2a2m 79 integer l_3a1m_sort,k_3a1m_sort 80 integer l_2a2m_aux,k_2a2m_aux 81 integer key_1a3m,offset_1a3m 82 integer l_1a3m,k_1a3m,size_1a3m 83 integer l_1a3m_sort,k_1a3m_sort 84 INTEGER IROW,ICOL,IRES 85 INTEGER INDEX_PAIR 86 integer l_4m,k_4m,size_4m,key_4m,offset_4m 87 integer sf3a1m_chunk,sf1a3m_chunk 88c 89 integer l_integral,l_coeff 90 integer k_integral,k_coeff 91 integer size_ic,size_icc,size_integral,size_coeff,max_na 92c 93 integer l_4af_offset,k_4af_offset,d_4af 94 integer k_2a2m_offset,l_2a2m_offset 95 integer sf_chunk,request 96 integer key_4af,offset_4af,size_4af 97 integer sf2a2m_chunk,key_2a2m,offset_2a2m 98 integer d_2a2m,size_2a2mf 99 character*255 filename 100c 101 integer size_max1,size_max2,ip_max1,ip_max2 102 integer chunk_max1,chunk_max2,i_from_x,i_to_x 103 integer isizef(5),imaxch(5) 104 integer d_max1,d_max2 105c 106 integer l_a4_offset,k_a4_offset 107 integer offset_bl4a,size_bl4a,a4_ini,a4_fin 108 integer offset_ismall,length_a4,key_ini 109 integer l_a3_offset,k_a3_offset 110 integer length_a3,offset_bl3a1m,size_bl3a1m 111 integer size_3a1mi,a3_ini,a3_fin 112c 113 integer l_a2_offset,k_a2_offset,length_a2 114 integer offset_bl2a2m,size_bl2a2m,a2_ini,a2_fin 115 integer size_1a3mi,size_2a2mi 116c 117 integer length_a1 118 integer l_a1_offset,k_a1_offset 119 integer offset_bl1a3m,size_bl1a3m,a1_ini,a1_fin 120c I/O improvements 121 integer handle_4a,handle_3a1m,handle_2a2m,handle_1a3m 122c 123c *** debug *** 124 double precision xxx,xmax 125c integer iqx 126 integer ierrcode1,ierrcode2 127c ************* 128c 129 logical parallel 130c 131 INTEGER length 132 INTEGER NXTASK 133 INTEGER next 134 INTEGER nprocs 135 INTEGER count 136 EXTERNAL NXTASK 137 integer nxtask 138 external nxtask 139 logical nodezero 140 logical idiskl 141c 142c 143c 144ccx idisk=.true. 145 if(idisk.eq.0) then 146 idiskl=.false. 147 else 148 idiskl=.true. 149 end if 150c 151 parallel = .true. 152c 153c 154 max_size_temp=imaxsize**4 155c 156c 157 do ii=1,1000 158 tot_zone(ii)=0.0d0 159 enddo 160 if(atpart.gt.1000) 161 & call errquit('tce_zones: atpart too big',1,MA_ERR) 162 sumx=0 163 do ii=1,atpart 164 tot_zone(ii)=sumx 165 sumx=sumx+nalength(ii) 166 enddo 167c 168 nodezero=(ga_nodeid().eq.0) 169c 170c *** debug *** 171c if(nodezero) then 172c write(6,*)'--------- NEW LOOP STRUCTURE --------' 173c call util_flush(6) 174c end if 175c ************* 176c 177c this module is called only if intorb = .true. 178c N is the number of correlated orbitals 179 N = nmo(1) - nfc(1) - nfv(1) 180 iha = nocc(1)-nfc(1) 181 ihb = nocc(ipol)-nfc(ipol) 182 ipa = nmo(1)-nocc(1)-nfv(1) 183 ipb = nmo(ipol)-nocc(ipol)-nfv(ipol) 184c 185 sf_chunk=(imaxsize+10)**4 186 sf3a1m_chunk=((imaxsize+10)**3)*tile_dim 187 sf2a2m_chunk=((tile_dim)**2)*((imaxsize+10)**2) 188 sf1a3m_chunk=((tile_dim)**3)*(imaxsize+10) 189c 190c l_integral and l_coeff local files are opened here 191c ATTENTION,ACHTUNG,UWAGA - "manually" defined 192cc size_ic=2*(imaxsize+10)**4 193cc size_ic=2*(imaxsize)**4 194 max_na=0 195 do ixi=1,atpart 196 if(nalength(ixi).gt.max_na) max_na=nalength(ixi) 197 enddo 198c it was 199cc size_ic=2*(max_na)**4+1 200c it is 201 size_ic=2*((max_na)**4)+1 202 size_icc=tile_dim*max_na 203c 204 if(nodezero) then 205 write(6,*)'2_EL_BATCH AND COEFF MATRIX' 206 write(6,*)'size_ic',size_ic 207 write(6,*)'size_icc',size_icc 208 write(6,*)'tile_dim',tile_dim 209 write(6,*)'max_na',max_na 210 call util_flush(6) 211 end if 212c 213 if (.not.ma_push_get(mt_dbl,size_ic,'l_int', 214 1 l_integral,k_integral)) 215 1 call errquit('tce_4s: MA problem l_int',0,MA_ERR) 216c 217 if (.not.ma_push_get(mt_dbl,size_icc,'l_coeff', 218 1 l_coeff,k_coeff)) 219 1 call errquit('tce_4s: MA problem l_coeff',0,MA_ERR) 220c 221 call tce_4a_offsetx2(l_4af_offset,k_4af_offset,size_4af) 222 call tce_3a1m_offsetx2(l_3a1m_offset,k_3a1m_offset,size_3a1mf) 223 call tce_2a2m_offsetx2(l_2a2m_offset,k_2a2m_offset,size_2a2mf) 224 call tce_1a3m_offsetx2(l_1a3m_offset,k_1a3m_offset,size_1a3mf) 225c 226c choose two largeast size_4af,size_3a1mf,size_2a2mf,size_1a3mf,size_2e 227c size_max1 > size_max2 228 size_max1=0 229 size_max2=0 230c 231 chunk_max1=0 232 chunk_max2=0 233c 234 isizef(1)=size_4af 235 isizef(2)=size_3a1mf 236 isizef(3)=size_2a2mf 237 isizef(4)=size_1a3mf 238 isizef(5)=size_2e 239c 240 imaxch(1)=sf_chunk 241 imaxch(2)=sf3a1m_chunk 242 imaxch(3)=sf2a2m_chunk 243 imaxch(4)=sf1a3m_chunk 244 imaxch(5)=tile_dim**4 245c 246 do i=1,5 247 if(isizef(i).gt.size_max1) then 248 size_max1=isizef(i) 249 ip_max1=i 250 chunk_max1=imaxch(i) 251 end if 252 enddo 253 do i=1,5 254 if(i.eq.ip_max1) go to 771 255 if(isizef(i).gt.size_max2) then 256 size_max2=isizef(i) 257 ip_max2=i 258 chunk_max2=imaxch(i) 259 end if 260 771 continue 261 enddo 262c 263c CHUNK_SIZE IS PROPERLY DEFINED HERE 264c 265 chunk_max1=size_ic 266 chunk_max2=size_ic 267c 268c 269c 270 if(nodezero) then 271 write(6,*)'INTERMEDIATE FILES' 272 write(6,*)'size_4af',size_4af 273 write(6,*)'size_3a1mf',size_3a1mf 274 write(6,*)'size_2a2mf',size_2a2mf 275 write(6,*)'size_1a3mf',size_1a3mf 276 write(6,*)'size_2e',size_2e 277 end if 278c 279c DISK USAGE ONLY !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 280c 281 if(idiskl) then 282 if(.not.parallel) 283 1 call errquit('sf only for parallel runs',1,DISK_ERR) 284cx if(parallel) call ga_sync() 285cx call util_file_name('max1',.false.,.false.,filename) 286cx if (sf_create(filename,dble(bytes)*dble(size_max1), 287cx 1 dble(bytes)*dble(size_max1),chunk_max1,d_max1).ne.0) 288cx 2 call errquit('4-index: sf problem a',0,DISK_ERR) 289c 290cx call util_file_name('max2',.false.,.false.,filename) 291cx if (sf_create(filename,dble(bytes)*dble(size_max2), 292cx 1 dble(bytes)*dble(size_max2),chunk_max2,d_max2).ne.0) 293cx 2 call errquit('4-index: sf problem b',0,DISK_ERR) 294c 295c call createfile(filename,d_v2,size_2e) 296c call ga_zero(d_v2) 297c 298 else 299 if(.not.parallel) 300 1 call errquit('sf only for parallel runs',1,DISK_ERR) 301 if(parallel) call ga_sync() 302 if(size_4af.le.size_max2) then 303 call createfile(filename,d_max1,size_max1) 304 call reconcilefile(d_max1,size_max1) 305 call ga_zero(d_max1) 306 call createfile(filename,d_max2,size_max2) 307 call reconcilefile(d_max2,size_max2) 308 call ga_zero(d_max2) 309 else 310 call createfile(filename,d_max2,size_max2) 311 call reconcilefile(d_max2,size_max2) 312 call ga_zero(d_max2) 313 call createfile(filename,d_max1,size_max1) 314 call reconcilefile(d_max1,size_max1) 315 call ga_zero(d_max1) 316 end if 317 end if 318c 319c 320c alpha orbitals only 321c 322 if (.not.ma_push_get(mt_dbl,nbf*(iha+ipa) 323 1 ,"sorted MO coeffs", 324 2 l_movecs_orb,k_movecs_orb)) 325 3 call errquit('tce_mo2e_zone: MA problem 1',0, 326 2 BASIS_ERR) 327 call dfill(nbf*(iha+ipa),0.0d0, dbl_mb(k_movecs_orb), 1) 328 do i=1,iha 329 do isigma1=1,nbf 330 dbl_mb(k_movecs_orb+(i-1)*nbf+isigma1-1)= 331 & dbl_mb(k_movecs_sorted+(i-1)*nbf+isigma1-1) 332 enddo 333 enddo 334 do i=iha+1,iha+ipa 335 do isigma1=1,nbf 336 dbl_mb(k_movecs_orb+(i-1)*nbf+isigma1-1)= 337 & dbl_mb(k_movecs_sorted+(i+ihb-1)*nbf+isigma1-1) 338 enddo 339 enddo 340c 341c 342 call int_mem_2e4c(work1,work2) 343 if (.not.ma_push_get(mt_dbl,work1,'work1',l_work1,k_work1)) 344 1 call errquit('tce_ao2e: MA problem work1',0,MA_ERR) 345 if (.not.ma_push_get(mt_dbl,work2,'work2',l_work2,k_work2)) 346 1 call errquit('tce_ao2e: MA problem work2',1,MA_ERR) 347c 348c 349c 350c 4af file formed here 351c 352ccx if(size_4af.le.size_max2) then 353ccx i_to_x=d_max2 354ccx else 355ccx i_to_x=d_max1 356ccx end if 357c 358c DISK CHANGES =============== 4A 359c opening handle_4a file 360c ** debug ** 361 if(nodezero) then 362 write(6,*)'before sf_create' 363 call util_flush(6) 364 end if 365c 366 if(idiskl) then 367 call util_file_name('4af',.false.,.false.,filename) 368 if (sf_create(filename,dble(bytes)*dble(size_4af), 369 1 dble(bytes)*dble(size_4af),chunk_max1,handle_4a).ne.0) 370 2 call errquit('4-index: sf problem opening 4af',0,DISK_ERR) 371 i_to_x=handle_4a 372 end if 373c 374 if(nodezero) then 375 write(6,*)'after sf_create' 376 call util_flush(6) 377 end if 378c ============================ 379 380 cpu1 = - util_cpusec() 381 wall1 = - util_wallsec() 382c 383 nprocs = GA_NNODES() 384 count = 0 385 next = NXTASK(nprocs, 1) 386 DO azone1 = 1,atpart !nu 387 DO azone2 = azone1,atpart !mu 388 DO azone3 = 1,atpart !sigma 389 DO azone4 = azone3,atpart !rho 390 IF (next.eq.count) THEN 391c --------------------------- 392 size_4a = nalength(azone1)*nalength(azone2)* 393 1 nalength(azone3)*nalength(azone4) 394 if(.not.ma_push_get(mt_dbl,size_4a,'4a',l_4a,k_4a)) 395 1 call errquit('tce_4af_zones1: MA problem',0,MA_ERR) 396 call dfill(size_4a, 0.0d0, dbl_mb(k_4a), 1) 397 shift_mu = 0 398 do mu = a2length(azone2)+1,a2length(azone2+1) 399 if (.not.bas_cn2bfr(ao_bas_han,mu,mu_lo,mu_hi)) 400 1 call errquit('tce_ao2e: basis fn range problem 1',0, 401 2 BASIS_ERR) 402 mu_range = mu_hi - mu_lo + 1 403 shift_nu = 0 404 do nu = a2length(azone1)+1,a2length(azone1+1) 405 if (.not.bas_cn2bfr(ao_bas_han,nu,nu_lo,nu_hi)) 406 1 call errquit('tce_ao2e: basis fn range problem 1',0, 407 2 BASIS_ERR) 408 nu_range = nu_hi - nu_lo + 1 409 shift_rho = 0 410 do rho = a2length(azone4)+1,a2length(azone4+1) 411 if (.not.bas_cn2bfr(ao_bas_han,rho,rho_lo,rho_hi)) 412 1 call errquit('tce_ao2e: basis fn range problem 1',0, 413 2 BASIS_ERR) 414 rho_range = rho_hi - rho_lo + 1 415 shift_sigma = 0 416 do sigma = a2length(azone3)+1,a2length(azone3+1) 417 if (.not.bas_cn2bfr(ao_bas_han,sigma,sigma_lo,sigma_hi)) 418 1 call errquit('tce_ao2e: basis fn range problem 1',0, 419 2 BASIS_ERR) 420 sigma_range = sigma_hi - sigma_lo + 1 421 if (schwarz_shell(rho,sigma)*schwarz_shell(mu,nu) 422 1 .ge. tol2e) then 423 call int_2e4c(ao_bas_han,mu,nu,ao_bas_han,rho,sigma, 424 1 work2,dbl_mb(k_work2),work1,dbl_mb(k_work1)) 425 i=0 426 do mu1 = 1,mu_range 427 do nu1 = 1,nu_range 428 do rho1 = 1,rho_range 429 do sigma1 = 1,sigma_range 430 i=i+1 431 inu1=nu1+shift_nu 432 isigma1=sigma1+shift_sigma 433 imu1=mu1+shift_mu 434 irho1=rho1+shift_rho 435c (isigma1,irho1|inu1, imu1) 436c differnce with other versiosn 437c here 438 ipos1=(((inu1-1)*nalength(azone2)+imu1-1)* 439 1 nalength(azone3)+isigma1-1)*nalength(azone4) 440 2 +irho1 441 dbl_mb(k_4a+ipos1-1)=dbl_mb(k_work1+i-1) 442 enddo 443 enddo 444 enddo 445 enddo 446 end if !schwarz screening 447 shift_sigma = shift_sigma + sigma_range 448 enddo !sigma 449 shift_rho = shift_rho + rho_range 450 enddo !rho 451 shift_nu = shift_nu + nu_range 452 enddo !nu 453 shift_mu = shift_mu + mu_range 454 enddo !mu 455c 456c fixing offsets and sf_writing 457 key_4af=azone4 - 1 + atpart * (azone3 - 1 + 458 & atpart * (azone2 - 1 + atpart * (azone1 - 1))) 459 call tce_hash(int_mb(k_4af_offset),key_4af,offset_4af) 460 if(idiskl) then !---- 461 ierrcode1=sf_write(i_to_x,dble(bytes)*dble(offset_4af), 462 1 dble(bytes)*dble(size_4a),dbl_mb(k_4a),request) 463ccx if (sf_write(i_to_x,dble(bytes)*dble(offset_4af), 464ccx 1 dble(bytes)*dble(size_4a),dbl_mb(k_4a),request).ne.0) 465ccx 1 THEN 466 if(ierrcode1.ne.0) then 467 write(6,8154)ierrcode1 468 write(6,8155)ga_nodeid(),key_4af,offset_4af,size_4a 469 call util_flush(6) 470 call errquit('zones put: sf problem2x1',1,DISK_ERR) 471 end if 472 ierrcode2=sf_wait(request) 473ccx if (sf_wait(request).ne.0) 474 if(ierrcode2.ne.0) then 475 write(6,8153)ierrcode2 476 call util_flush(6) 477 call errquit('zones put: sf problem3x1',2,DISK_ERR) 478 end if 479 else 480 call ga_put(i_to_x,offset_4af+1,offset_4af+size_4a,1,1, 481 1 dbl_mb(k_4a),1) 482 end if !---- 483c closing l_4a file 484 if (.not.ma_pop_stack(l_4a)) 485 1 call errquit('tce_mo2e_4af2: l_4a',15,MA_ERR) 486c --------------------------- 487 next = NXTASK(nprocs, 1) 488 END IF 489 count = count + 1 490 ENDDO !azone4 491 ENDDO !azone3 492 ENDDO !azone2 493 ENDDO !azone1 494c 495c 496 call ga_sync() 497 next = nxtask(-nprocs,1) 498c 499c 500c 501 if(nodezero) then 502 write(6,*)'STEP2 4index' 503 call util_flush(6) 504 end if 505c 506c 507c 508c 509c 510ccx i_from_x=i_to_x 511ccx if(i_from_x.eq.d_max1) then 512ccx i_to_x=d_max2 513c size check-out 514ccx if(size_max2.lt.size_3a1mf) 515ccx 1 call errquit('tce_4s: 1a',15,MA_ERR) 516ccx else 517ccx i_to_x=d_max1 518c size check-out 519ccx if(size_max1.lt.size_3a1mf) 520ccx 1 call errquit('tce_4s: 1b',15,MA_ERR) 521ccx end if 522c 523c DISK CHANGES =============== 3A1M 524c closing/reopening handle_4a file 525c opening handle_3a1m file 526 if(idiskl) then 527c 528 if(sf_rwtor(handle_4a).ne.0) 529 2 call errquit('4-index: sf_rwtor 4a',0,DISK_ERR) 530c *** debug *** 531c if(nodezero) then 532c write(6,*)'STEP2 tests' 533c write(6,*)'handle_4a',handle_4a 534c call util_flush(6) 535c endif 536c ************* 537c 538 call util_file_name('3a1mf',.false.,.false.,filename) 539 if (sf_create(filename,dble(bytes)*dble(size_3a1mf), 540 1 dble(bytes)*dble(size_3a1mf),chunk_max1,handle_3a1m).ne.0) 541 2 call errquit('4-index: sf problem opening 3a1mf',0,DISK_ERR) 542 i_from_x=handle_4a 543 i_to_x=handle_3a1m 544c *** debug *** 545c if(nodezero) then 546c write(6,*)'i_from_x=',i_from_x 547c write(6,*)'i_to_x=',i_to_x 548c write(6,*)'handle_4a=',handle_4a 549c write(6,*)'handle_3a1m=',handle_3a1m 550c write(6,*)'size_3a1mf=',size_3a1mf 551c write(6,*)'size_4af=',size_4af 552c call util_flush(6) 553c end if 554c ************* 555 end if 556c ============================ 557 558c 559c 3A1M part here 560c 561 nprocs = GA_NNODES() 562 count = 0 563 next = NXTASK(nprocs, 1) 564c 565 DO azone1 = 1,atpart !nu 566 DO azone2 = azone1,atpart !mu 567 DO g4b = 1,noa+nva !k 568 DO azone3 = 1,atpart !sigma 569ccx DO g4b = 1,noa+nva !k 570c 571 IF (next.eq.count) THEN 572c 573 size_3a1m=int_mb(k_range_alpha+g4b-1) 574 1 *nalength(azone1)*nalength(azone2)*nalength(azone3) 575 if (.not.ma_push_get(mt_dbl,size_3a1m,'3a1m',l_3a1m,k_3a1m)) 576 1 call errquit('tce_4ind: step1_1',0,MA_ERR) 577 call dfill(size_3a1m, 0.0d0, dbl_mb(k_3a1m), 1) 578c 579c routine for offset 580c 581 call a4_offset(azone1,azone2,azone3,size_ic, 582 & l_a4_offset,k_a4_offset) 583c 584 length_a4=int_mb(k_a4_offset) 585 DO ii=1,length_a4 586 a4_ini = int_mb(k_a4_offset+ii) 587 a4_fin = int_mb(k_a4_offset+length_a4+ii) 588 key_ini=a4_ini - 1 + atpart * (azone3 - 1 + 589 & atpart * (azone2 - 1 + atpart * (azone1 - 1))) 590 call tce_hash(int_mb(k_4af_offset),key_ini,offset_bl4a) 591c 592 size_bl4a = int_mb(k_a4_offset+2*length_a4+ii) 593 if(idiskl) then !--------- 594 ierrcode1=sf_read(i_from_x,dble(bytes)*dble(offset_bl4a), 595 1 dble(bytes)*dble(size_bl4a),dbl_mb(k_integral), 596 1 request) 597ccx if (sf_read(i_from_x,dble(bytes)*dble(offset_bl4a), 598ccx 1 dble(bytes)*dble(size_bl4a),dbl_mb(k_integral), 599ccx 1 request).ne.0) THEN 600 if(ierrcode1.ne.0) THEN 601 write(6,*)'STEP2_step1' 602 write(6,8153)ierrcode1 603 write(6,8155)ga_nodeid(),key_ini,offset_bl4a,size_bl4a 604 write(6,*)'length_a4',length_a4 605 write(6,*)'ii=',ii 606 write(6,*)'i_from_x',i_from_x 607 call util_flush(6) 608 call errquit('zones put: sf problem2x2',1,DISK_ERR) 609 end if 610 ierrcode2=sf_wait(request) 611ccx if (sf_wait(request).ne.0) 612 if (ierrcode2.ne.0) then 613 write(6,8154)ierrcode2 614 call util_flush(6) 615 call errquit('zones put: sf problem3x2',2,DISK_ERR) 616 end if 617 else 618 call ga_get(i_from_x,offset_bl4a+1,offset_bl4a+size_bl4a,1,1, 619 1 dbl_mb(k_integral),1) 620 end if !-------------- 621c 622 offset_ismall=0 623 do i=a4_ini,a4_fin ! over small i ------------------ 624 tot_azone4_sh=tot_zone(i) 625 j=0 626 do irho1 = 1,nalength(i) 627 do igi4=1,int_mb(k_range_alpha+g4b-1) 628 j=j+1 629 ipos1=(int_mb(k_offset_alpha+g4b-1)+igi4-1)*nbf+tot_azone4_sh 630 & +irho1 631 dbl_mb(k_coeff+j-1)=dbl_mb(k_movecs_orb+ipos1-1) 632 enddo 633 enddo 634c C(g4b irho)* v(irho,isigma,imu,inu) 635 call dgemm('N','N', 636 & int_mb(k_range_alpha+g4b-1), 637 & nalength(azone1)*nalength(azone2)*nalength(azone3), 638 & nalength(i), 639 & 1.0d0,dbl_mb(k_coeff),int_mb(k_range_alpha+g4b-1), 640 & dbl_mb(k_integral+offset_ismall),nalength(i),1.0d0, 641 & dbl_mb(k_3a1m),int_mb(k_range_alpha+g4b-1)) 642c 643 offset_ismall=offset_ismall+nalength(azone1)*nalength(azone2)* 644 1 nalength(azone3)*nalength(i) 645 enddo ! over small i ---------------------- 646c 647 ENDDO ! over ii 648 if (.not.ma_pop_stack(l_a4_offset)) 649 1 call errquit('tce_mo2e_trans_zones: uu l_2g2z',15,MA_ERR) 650c 651c 652c 653c 654 DO azone4 = 1, azone3-1 ! second part ---- 655c getting piece of 4a 656 size_4a = nalength(azone1)*nalength(azone2)* 657 1 nalength(azone3)*nalength(azone4) 658 if(azone3.le.azone4) then 659 key_4af=azone4 - 1 + atpart * (azone3 - 1 + 660 & atpart * (azone2 - 1 + atpart * (azone1 - 1))) 661 else 662 key_4af=azone3 - 1 + atpart * (azone4 - 1 + 663 & atpart * (azone2 - 1 + atpart * (azone1 - 1))) 664 end if 665 call tce_hash(int_mb(k_4af_offset),key_4af,offset_4af) 666 if(idiskl) then !--------- 667 ierrcode1=sf_read(i_from_x,dble(bytes)*dble(offset_4af), 668 1 dble(bytes)*dble(size_4a),dbl_mb(k_integral),request) 669ccx if (sf_read(i_from_x,dble(bytes)*dble(offset_4af), 670ccx 1 dble(bytes)*dble(size_4a),dbl_mb(k_integral),request).ne.0) 671ccx 1 THEN 672 if(ierrcode1.ne.0) then 673 write(6,*)'STEP2-step1b' 674 write(6,8153)ierrcode1 675 write(6,8155)ga_nodeid(),key_4af,offset_4af,size_4a 676 call errquit('zones put: sf problem2x3',1,DISK_ERR) 677 end if 678 ierrcode2=sf_wait(request) 679ccx if (sf_wait(request).ne.0) 680 if(ierrcode2.ne.0) then 681 write(6,8154)ierrcode2 682 call util_flush(6) 683 call errquit('zones put: sf problem3x3',2,DISK_ERR) 684 end if 685 else 686 call ga_get(i_from_x,offset_4af+1,offset_4af+size_4a,1,1, 687 1 dbl_mb(k_integral),1) 688 end if !-------------- 689c 690 if(azone4.lt.azone3) then 691 if(.not.ma_push_get(mt_dbl,size_4a,'4a_s',l_4a_sort,k_4a_sort)) 692 1 call errquit('tce_4af_zones1: MA problem',0,MA_ERR) 693 CALL TCE_SORT_4KG_(dbl_mb(k_integral),dbl_mb(k_4a_sort), 694 & nalength(azone3),nalength(azone4), 695 & nalength(azone2),nalength(azone1), 696 &2,1,3,4,1.0d0) 697 do i=1,size_4a 698 dbl_mb(k_integral+i-1)=dbl_mb(k_4a_sort+i-1) 699 enddo 700 if (.not.ma_pop_stack(l_4a_sort)) 701 1 call errquit('tce_mo2e_trans_zones: uu l_2g2z',15,MA_ERR) 702 end if 703c 704c C(g4b irho) ==> l_2g2z 705c 706 size_amc=int_mb(k_range_alpha+g4b-1) * 707 1 nalength(azone4) 708 i=0 709 tot_azone4_sh=tot_zone(azone4) 710 do irho1 = 1,nalength(azone4) 711 do igi4=1,int_mb(k_range_alpha+g4b-1) 712 i=i+1 713 ipos1=(int_mb(k_offset_alpha+g4b-1)+igi4-1)*nbf+tot_azone4_sh 714 & +irho1 715 dbl_mb(k_coeff+i-1)=dbl_mb(k_movecs_orb+ipos1-1) 716 enddo 717 enddo 718c C(g4b irho)* v(irho,isigma,imu,inu) 719 call dgemm('N','N', 720 & int_mb(k_range_alpha+g4b-1), 721 & nalength(azone1)*nalength(azone2)*nalength(azone3), 722 & nalength(azone4), 723 & 1.0d0,dbl_mb(k_coeff),int_mb(k_range_alpha+g4b-1), 724 & dbl_mb(k_integral),nalength(azone4),1.0d0, 725 & dbl_mb(k_3a1m),int_mb(k_range_alpha+g4b-1)) 726c 727 ENDDO 728c write to file 729 key_3a1m = azone3-1+atpart*(g4b-1+(noa+nva)* 730 & (azone2 - 1 + atpart * (azone1 - 1))) 731 call tce_hash(int_mb(k_3a1m_offset),key_3a1m,offset_3a1m) 732 if(idiskl) then ! -------- 733 ierrcode1=sf_write(i_to_x,dble(bytes)*dble(offset_3a1m), 734 r dble(bytes)*dble(size_3a1m),dbl_mb(k_3a1m),request) 735 if (ierrcode1.ne.0) then 736 write(6,*)'WRITE STEP3' 737 write(6,8156) ierrcode1 738 write(6,*)'offset_3a1m=',offset_3a1m 739 write(6,*)'size_3a1m=',size_3a1m 740 write(6,*)'i_to_x',i_to_x 741 call util_flush(6) 742 call errquit('zones put: sf problem21-b',1,DISK_ERR) 743 end if 744 ierrcode2=sf_wait(request) 745ccx if (sf_wait(request).ne.0) 746 if(ierrcode2.ne.0) then 747 write(6,8157)ierrcode2 748 call errquit('zones put: sf problem31-b',2,DISK_ERR) 749 end if 750 else 751 call ga_put(i_to_x,offset_3a1m+1,offset_3a1m+size_3a1m,1,1, 752 1 dbl_mb(k_3a1m),1) 753 end if ! --------- 754c 755 if (.not.ma_pop_stack(l_3a1m)) 756 1 call errquit('tce_mo2e_trans_zones: abc-MA problem',15,MA_ERR) 757c 758 next = NXTASK(nprocs, 1) 759 END IF 760 count = count + 1 761 ENDDO 762 ENDDO 763 ENDDO 764 ENDDO 765c 3A1P fully done here 766 call ga_sync() 767 next = nxtask(-nprocs,1) 768c 769c 770c 771c 772 if(nodezero) then 773 write(6,*)'STEP3 4index' 774 call util_flush(6) 775 end if 776c 777c 778c 779c 780c 781c 2A2M part here 782c 783ccx i_from_x=i_to_x 784ccx if(i_from_x.eq.d_max1) then 785ccx i_to_x=d_max2 786c size check-out 787ccx if(size_max2.lt.size_2a2mf) 788ccx 1 call errquit('tce_4s: 2a',15,MA_ERR) 789ccx else 790ccx i_to_x=d_max1 791c size check-out 792ccx if(size_max1.lt.size_2a2mf) 793ccx 1 call errquit('tce_4s: 2b',15,MA_ERR) 794ccx end if 795c 796c 797c DISK CHANGES =============== 2A2M 798c destroying handle_4a 799c closing/reopening handle_3a1m file 800c opening handle_2a2m file 801 if(idiskl) then 802c 803 if (sf_destroy(handle_4a).ne.0) 804 1 call errquit('tce_sf_destroy: sf problem handle_4a',15,MA_ERR) 805c 806 if(sf_rwtor(handle_3a1m).ne.0) 807 2 call errquit('4-index: sf_rwtor 3a1m',0,DISK_ERR) 808c *** debug *** 809c if(nodezero) then 810c write(6,*)'handle_3a1m',handle_3a1m 811c end if 812c ************* 813c 814 call util_file_name('2a2mf',.false.,.false.,filename) 815 if (sf_create(filename,dble(bytes)*dble(size_2a2mf), 816 1 dble(bytes)*dble(size_2a2mf),chunk_max1,handle_2a2m).ne.0) 817 2 call errquit('4-index: sf problem opening 2a2mf',0,DISK_ERR) 818 i_from_x=handle_3a1m 819 i_to_x=handle_2a2m 820 end if 821c *** debug *** 822 if(nodezero) then 823 write(6,*)'STEP3 test' 824c write(6,*)'handle_2a2m=',handle_2a2m 825c write(6,*)'size_2a2mf=',size_2a2mf 826c write(6,*)'i_from_x=',i_from_x 827c write(6,*)'i_to_x=',i_to_x 828 call util_flush(6) 829 end if 830c ************* 831c ============================ 832c 833 nprocs = GA_NNODES() 834 count = 0 835 next = NXTASK(nprocs, 1) 836c 837 DO g3b = 1,noa+nva !k 838 DO g4b = g3b,noa+nva !l 839 DO azone1 = 1,atpart !nu 840 DO azone2 = azone1,atpart !mu 841ccx DO g3b = 1,noa+nva !k 842ccx DO g4b = g3b,noa+nva !l 843c 844 IF (next.eq.count) THEN 845c 846 size_2a2m=int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 847 1 *nalength(azone1)*nalength(azone2) 848 if (.not.ma_push_get(mt_dbl,size_2a2m,'2a2m',l_2a2m,k_2a2m)) 849 1 call errquit('tce_4ind: step1_1',0,MA_ERR) 850 call dfill(size_2a2m, 0.0d0, dbl_mb(k_2a2m), 1) 851c 852 call a3_offset(azone1,azone2,g3b,size_ic, 853 & l_a3_offset,k_a3_offset) 854 length_a3=int_mb(k_a3_offset) 855c 856 DO ii=1,length_a3 ! ii---loop 857 858 a3_ini = int_mb(k_a3_offset+ii) 859 a3_fin = int_mb(k_a3_offset+length_a3+ii) 860 key_ini = a3_ini-1+atpart*(g3b-1+(noa+nva)* 861 & (azone2 - 1 + atpart * (azone1 - 1))) 862 call tce_hash(int_mb(k_3a1m_offset),key_ini,offset_bl3a1m) 863c 864 size_bl3a1m = int_mb(k_a3_offset+2*length_a3+ii) 865c 866c 867 if(idiskl) then ! ------- 868 ierrcode1=sf_read(i_from_x,dble(bytes)*dble(offset_bl3a1m), 869 2 dble(bytes)*dble(size_bl3a1m), 870 2 dbl_mb(k_integral),request) 871 if(ierrcode1.ne.0) then 872 write(6,*)'STEP3-step1a' 873 write(6,8153) ierrcode1 874 write(6,*)'offset_bl3a1m=',offset_bl3a1m 875 write(6,*)'size_bl3a1m=',size_bl3a1m 876 write(6,*)'ii=',ii 877 write(6,*)'length_a3=',length_a3 878 call util_flush(6) 879 call errquit('zones put: sf problem21-b',1,DISK_ERR) 880 end if 881 ierrcode2=sf_wait(request) 882 if(ierrcode2.ne.0) then 883 write(6,*)'STEP3-step1a' 884 write(6,8154) ierrcode2 885 call util_flush(6) 886 call errquit('zones put: sf problem31-b',2,DISK_ERR) 887 end if 888 else 889 call ga_get(i_from_x,offset_bl3a1m+1,offset_bl3a1m+size_bl3a1m, 890 1 1,1,dbl_mb(k_integral),1) 891 end if ! ------- 892 offset_ismall=0 893 do i=a3_ini,a3_fin ! over small i ------------------ 894 size_3a1mi=nalength(azone1)*nalength(azone2)* 895 1 nalength(i)*int_mb(k_range_alpha+g3b-1) 896 if (.not.ma_push_get(mt_dbl,size_3a1mi,'3a1m', 897 1 l_3a1m_sort,k_3a1m_sort)) 898 1 call errquit('tce_4ind: step1_1',0,MA_ERR) 899 CALL TCE_SORT_4KG_(dbl_mb(k_integral+offset_ismall), 900 & dbl_mb(k_3a1m_sort), 901 & int_mb(k_range_alpha+g3b-1),nalength(i), 902 & nalength(azone2),nalength(azone1), 903 &2,1,3,4,1.0d0) 904 do j=1,size_3a1mi 905 dbl_mb(k_integral+offset_ismall+j-1)=dbl_mb(k_3a1m_sort+j-1) 906 enddo 907 if (.not.ma_pop_stack(l_3a1m_sort)) 908 1 call errquit('tce_mo2e_trans_zones: abc-MA problem',15,MA_ERR) 909c C(g4b isigma) 910 j=0 911 tot_azone3_sh=tot_zone(i) 912 do isigma1 = 1,nalength(i) 913 do igi4=1,int_mb(k_range_alpha+g4b-1) 914 j=j+1 915 ipos1=(int_mb(k_offset_alpha+g4b-1)+igi4-1)*nbf+tot_azone3_sh 916 & +isigma1 917 dbl_mb(k_coeff+j-1)=dbl_mb(k_movecs_orb+ipos1-1) 918 enddo 919 enddo 920c C(g4b isigma)* v(isigma,g3b,imu,inu) 921 call dgemm('N','N', 922 & int_mb(k_range_alpha+g4b-1), 923 & nalength(azone1)*nalength(azone2)*int_mb(k_range_alpha+g3b-1), 924 & nalength(i), 925 & 1.0d0,dbl_mb(k_coeff),int_mb(k_range_alpha+g4b-1), 926 & dbl_mb(k_integral+offset_ismall),nalength(i),1.0d0, 927 & dbl_mb(k_2a2m),int_mb(k_range_alpha+g4b-1)) 928c 929 offset_ismall=offset_ismall+nalength(azone1)*nalength(azone2)* 930 1 nalength(i)*int_mb(k_range_alpha+g3b-1) 931 enddo ! over small i ---------------------- 932 ENDDO ! ii---loop 933c 934 if (.not.ma_pop_stack(l_a3_offset)) 935 1 call errquit('tce_mo2e_trans_zones: uu l_2g2z',15,MA_ERR) 936c write to file 937 key_2a2m= azone2 - 1 + atpart*( azone1 - 1 + 938 & atpart*(g4b - 1 + (noa+nva) * (g3b-1))) 939 call tce_hash(int_mb(k_2a2m_offset),key_2a2m,offset_2a2m) 940 if(idiskl) then ! ------- 941 ierrcode1=sf_write(i_to_x,dble(bytes)*dble(offset_2a2m), 942 2 dble(bytes)*dble(size_2a2m),dbl_mb(k_2a2m),request) 943 if(ierrcode1.ne.0) then 944 write(6,*)'STEP3 write' 945 write(6,8156) ierrcode1 946 write(6,*)'offset_2a2m=',offset_2a2m 947 write(6,*)'size_2a2m=',size_2a2m 948 write(6,*)'key_2a2m=',key_2a2m 949 write(6,*)'i_to_x=',i_to_x 950 call util_flush(6) 951 call errquit('zones put: sf problem21-b',1,DISK_ERR) 952 end if 953 ierrcode2=sf_wait(request) 954 if(ierrcode2.ne.0) then 955 call errquit('zones put: sf problem31-b',2,DISK_ERR) 956 end if 957 else 958 call ga_put(i_to_x,offset_2a2m+1,offset_2a2m+size_2a2m,1,1, 959 1 dbl_mb(k_2a2m),1) 960 end if ! -------- 961c 962 if (.not.ma_pop_stack(l_2a2m)) 963 1 call errquit('tce_mo2e_trans_zones: abc-MA problem',15,MA_ERR) 964c 965 next = NXTASK(nprocs, 1) 966 END IF 967 count = count + 1 968 ENDDO 969 ENDDO 970 ENDDO 971 ENDDO 972c 2A2P fully done here 973 call ga_sync() 974 next = nxtask(-nprocs,1) 975c 976c 977c 978c 979c 980 if(nodezero) then 981 write(6,*)'STEP4 4index' 982 call util_flush(6) 983 end if 984c 985c 986c 987c 988c 989c 990c 1A3M part here ([g4b][g3b]|[mu][nu]) =>([g4b][g3b]|[g2b][nu]) 991c [mu]>[nu] (azone2>=azone1) 992c 993c 994cx i_from_x=i_to_x 995cx if(i_from_x.eq.d_max1) then 996cx i_to_x=d_max2 997c size check-out 998cx if(size_max2.lt.size_1a3mf) 999cx 1 call errquit('tce_4s: 3a',15,MA_ERR) 1000cx else 1001cx i_to_x=d_max1 1002c size check-out 1003cx if(size_max1.lt.size_1a3mf) 1004cx 1 call errquit('tce_4s: 3b',15,MA_ERR) 1005cx end if 1006c 1007c 1008c 1009c 1010c DISK CHANGES =============== 1A3M 1011c destroying handle_3a1m 1012c closing/reopening handle_2a2m file 1013c opening handle_1a3m file 1014 if(idiskl) then 1015c 1016 if (sf_destroy(handle_3a1m).ne.0) 1017 1 call errquit('tce_sf_destroy: sf problem handle_3a1m',15,MA_ERR) 1018c 1019 if(sf_rwtor(handle_2a2m).ne.0) 1020 2 call errquit('4-index: sf_rwtor 2a2m',0,DISK_ERR) 1021c 1022 call util_file_name('1a3mf',.false.,.false.,filename) 1023 if (sf_create(filename,dble(bytes)*dble(size_1a3mf), 1024 1 dble(bytes)*dble(size_1a3mf),chunk_max1,handle_1a3m).ne.0) 1025 2 call errquit('4-index: sf problem opening 1a3mf',0,DISK_ERR) 1026 i_from_x=handle_2a2m 1027 i_to_x=handle_1a3m 1028 end if 1029c ============================ 1030c 1031c 1032c 1033c 1034 nprocs = GA_NNODES() 1035 count = 0 1036 next = NXTASK(nprocs, 1) 1037c 1038 DO g3b = 1,noa+nva !k 1039 DO g4b = g3b,noa+nva !l 1040 DO g2b = 1,noa+nva ! 1041 DO azone1 = 1,atpart !nu 1042ccx DO g3b = 1,noa+nva !k 1043ccx DO g4b = g3b,noa+nva !l 1044c 1045 IF (next.eq.count) THEN 1046 size_1a3m=int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1047 1 *int_mb(k_range_alpha+g2b-1)*nalength(azone1) 1048 if (.not.ma_push_get(mt_dbl,size_1a3m,'1a3m',l_1a3m,k_1a3m)) 1049 1 call errquit('tce_4ind: step1_1',0,MA_ERR) 1050 call dfill(size_1a3m, 0.0d0, dbl_mb(k_1a3m), 1) 1051c azone2 >= azone1 - PART1 1052c 1053c 1054 call a2_offset(azone1,g4b,g3b,size_ic, 1055 & l_a2_offset,k_a2_offset) 1056c 1057 length_a2=int_mb(k_a2_offset) 1058 DO ii=1,length_a2 !========================================= 1059 a2_ini = int_mb(k_a2_offset+ii) 1060 a2_fin = int_mb(k_a2_offset+length_a2+ii) 1061 key_ini=a2_ini - 1 + atpart*( azone1 - 1 + 1062 & atpart*(g4b - 1 + (noa+nva) * (g3b-1))) 1063 call tce_hash(int_mb(k_2a2m_offset),key_ini,offset_bl2a2m) 1064c 1065 size_bl2a2m = int_mb(k_a2_offset+2*length_a2+ii) 1066 if(idiskl) then !--------- 1067 if (sf_read(i_from_x,dble(bytes)*dble(offset_bl2a2m), 1068 1 dble(bytes)*dble(size_bl2a2m),dbl_mb(k_integral), 1069 1 request).ne.0) THEN 1070 write(6,8155)ga_nodeid(),key_4af,offset_4af,size_4a 1071 call errquit('zones put: sf problem2x4',1,DISK_ERR) 1072 end if 1073 if (sf_wait(request).ne.0) 1074 1 call errquit('zones put: sf problem3x4',2,DISK_ERR) 1075 else 1076 call ga_get(i_from_x,offset_bl2a2m+1, 1077 1 offset_bl2a2m+size_bl2a2m,1,1,dbl_mb(k_integral),1) 1078 end if !-------------- 1079c 1080c *** debug *** 1081c xmax=0.0d0 1082c do i=1,size_bl2a2m 1083c xmax=xmax+dbl_mb(k_integral+i-1) 1084c enddo 1085c if (util_print('mo2e',print_debug)) write(6,411) 1086c & azone1,g2b,g4b,g3b,a2_ini,a2_fin,size_bl2a2m,key_ini,xmax 1087c 411 format(' 1A3M part1 ',4i3,2x,2i4,2x,i9,2x,i9,f17.8) 1088c call util_flush(6) 1089c ************* 1090c 1091 offset_ismall=0 1092 do i=a2_ini,a2_fin ! over small i ------------------ 1093ccc TRANSPOZYCJA 1094 size_2a2mi=nalength(azone1)*nalength(i)* 1095 1 int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1096 if (.not.ma_push_get(mt_dbl,size_2a2mi,'2a2m_sort', 1097 1 l_2a2m_aux,k_2a2m_aux)) 1098 1 call errquit('tce_4ind: step1_1',0,MA_ERR) 1099 CALL TCE_SORT_4KG_(dbl_mb(k_integral+offset_ismall), 1100 & dbl_mb(k_2a2m_aux), 1101 & int_mb(k_range_alpha+g4b-1),int_mb(k_range_alpha+g3b-1), 1102 & nalength(i),nalength(azone1), 1103 & 1,2,4,3,1.0d0) 1104 do j=1,size_2a2mi 1105 dbl_mb(k_integral+offset_ismall+j-1)=dbl_mb(k_2a2m_aux+j-1) 1106 enddo 1107 if (.not.ma_pop_stack(l_2a2m_aux)) 1108 1 call errquit('tce_mo2e_trans_zones: abc-MA problem',15,MA_ERR) 1109c C(azone2 g2b) 1110 j=0 1111 tot_azone2_sh=tot_zone(i) 1112 do igi2=1,int_mb(k_range_alpha+g2b-1) 1113 do imu1 = 1,nalength(i) 1114 j=j+1 1115 ipos1=(int_mb(k_offset_alpha+g2b-1)+igi2-1)*nbf+tot_azone2_sh 1116 & +imu1 1117 dbl_mb(k_coeff+j-1)=dbl_mb(k_movecs_orb+ipos1-1) 1118 enddo 1119 enddo 1120c v(g4b,g3b,inu,imu) C(imu,g2b) 1121 size_1a3mi=int_mb(k_range_alpha+g4b-1) 1122 & *int_mb(k_range_alpha+g3b-1)*nalength(azone1) 1123 & *int_mb(k_range_alpha+g2b-1) 1124 if (.not.ma_push_get(mt_dbl,size_1a3mi,'1a3m', 1125 1 l_1a3m_sort,k_1a3m_sort)) 1126 1 call errquit('tce_4ind: step1_1',0,MA_ERR) 1127c dgemm 1128 call dgemm('N','N', 1129 & int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1130 & *nalength(azone1), 1131 & int_mb(k_range_alpha+g2b-1), 1132 & nalength(i), 1133 & 1.0d0,dbl_mb(k_integral+offset_ismall), 1134 & int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1135 & *nalength(azone1), 1136 & dbl_mb(k_coeff),nalength(i),0.0d0, 1137 & dbl_mb(k_1a3m_sort), 1138 & int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1139 & *nalength(azone1)) 1140c transpose 1141 j=0 1142 do inu1=1,nalength(azone1) 1143 do igi2=1,int_mb(k_range_alpha+g2b-1) 1144 do igi3=1,int_mb(k_range_alpha+g3b-1) 1145 do igi4=1,int_mb(k_range_alpha+g4b-1) 1146 j=j+1 1147 ipos1=igi4+int_mb(k_range_alpha+g4b-1)*( 1148 & igi3-1+int_mb(k_range_alpha+g3b-1)*( 1149 & inu1-1+nalength(azone1)*(igi2-1))) 1150 dbl_mb(k_1a3m+j-1)=dbl_mb(k_1a3m+j-1)+ 1151 & dbl_mb(k_1a3m_sort+ipos1-1) 1152 enddo 1153 enddo 1154 enddo 1155 enddo 1156 if (.not.ma_pop_stack(l_1a3m_sort)) 1157 1 call errquit('tce_mo2e_trans_zones: abc-MA problem',15,MA_ERR) 1158 offset_ismall=offset_ismall+nalength(azone1)*nalength(i)* 1159 1 int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1160 enddo ! over small i ---------------------- 1161 ENDDO ! over ii !=========================================== 1162 if (.not.ma_pop_stack(l_a2_offset)) 1163 1 call errquit('tce_mo2e_trans_zones: uu l_2g2z',15,MA_ERR) 1164c 1165c 1166c azone2 < azone1 - PART2 1167 1168 DO azone2 = 1,azone1-1 ! PART2 two cases here only one is followed 1169 size_2a2m=int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1170 1 *nalength(azone1)*nalength(azone2) 1171 if (.not.ma_push_get(mt_dbl,size_2a2m,'2a2m',l_2a2m,k_2a2m)) 1172 1 call errquit('tce_4ind: step1_1',0,MA_ERR) 1173 if(azone2.ge.azone1) then 1174 key_2a2m= azone2 - 1 + atpart*( azone1 - 1 + 1175 & atpart*(g4b - 1 + (noa+nva) * (g3b-1))) 1176 call tce_hash(int_mb(k_2a2m_offset),key_2a2m,offset_2a2m) 1177 if(idiskl) then !------ 1178 if (sf_read(i_from_x,dble(bytes)*dble(offset_2a2m), 1179 r dble(bytes)*dble(size_2a2m),dbl_mb(k_2a2m),request).ne.0) 1180 2 call errquit('zones put: sf problem21-b',1,DISK_ERR) 1181 if (sf_wait(request).ne.0) 1182 1 call errquit('zones put: sf problem31-b',2,DISK_ERR) 1183 else 1184 call ga_get(i_from_x,offset_2a2m+1,offset_2a2m+size_2a2m,1,1, 1185 1 dbl_mb(k_2a2m),1) 1186 end if !----- 1187 else 1188 key_2a2m= azone1 - 1 + atpart*( azone2 - 1 + 1189 & atpart*(g4b - 1 + (noa+nva) * (g3b-1))) 1190 call tce_hash(int_mb(k_2a2m_offset),key_2a2m,offset_2a2m) 1191 if(idiskl) then ! ------ 1192 if (sf_read(i_from_x,dble(bytes)*dble(offset_2a2m), 1193 r dble(bytes)*dble(size_2a2m),dbl_mb(k_2a2m),request).ne.0) 1194 2 call errquit('zones put: sf problem21-b',1,DISK_ERR) 1195 if (sf_wait(request).ne.0) 1196 1 call errquit('zones put: sf problem31-b',2,DISK_ERR) 1197c *** debug *** 1198c xmax=0.0d0 1199c do i=1,size_2a2m 1200c xmax=xmax+dbl_mb(k_2a2m+i-1) 1201c enddo 1202c if (util_print('mo2e',print_debug)) write(6,412) 1203c & azone1,g2b,g4b,g3b,azone2,size_2a2m,key_2a2m,xmax 1204c 412 format(' 1A3M part2 ',4i3,2x,i4,2x,i9,2x,i9,2x,f17.8) 1205c call util_flush(6) 1206c ************* 1207 else 1208 call ga_get(i_from_x,offset_2a2m+1,offset_2a2m+size_2a2m,1,1, 1209 1 dbl_mb(k_2a2m),1) 1210 end if ! ------ 1211 end if ! xxxxxx 1212 if(azone2.ge.azone1) then 1213ccc TRANSPOZYCJA 1214 if (.not.ma_push_get(mt_dbl,size_2a2m,'2a2m_sort', 1215 1 l_2a2m_aux,k_2a2m_aux)) 1216 1 call errquit('tce_4ind: step1_1',0,MA_ERR) 1217 CALL TCE_SORT_4KG_(dbl_mb(k_2a2m),dbl_mb(k_2a2m_aux), 1218 & int_mb(k_range_alpha+g4b-1),int_mb(k_range_alpha+g3b-1), 1219 & nalength(azone2),nalength(azone1), 1220 & 1,2,4,3,1.0d0) 1221 do i=1,size_2a2m 1222 dbl_mb(k_2a2m+i-1)=dbl_mb(k_2a2m_aux+i-1) 1223 enddo 1224 if (.not.ma_pop_stack(l_2a2m_aux)) 1225 1 call errquit('tce_mo2e_trans_zones: abc-MA problem',15,MA_ERR) 1226 end if 1227c 1228c C(imu g2b) ==> l_2g2z 1229c 1230 size_amc=int_mb(k_range_alpha+g2b-1) * 1231 1 nalength(azone2) 1232 if (.not.ma_push_get(mt_dbl,size_amc,'2g2z',l_2g2z,k_2g2z)) 1233 1 call errquit('tce_r2_divide3: xx-MA problem',0,MA_ERR) 1234 i=0 1235 tot_azone2_sh=tot_zone(azone2) 1236 do igi2=1,int_mb(k_range_alpha+g2b-1) 1237 do imu1 = 1,nalength(azone2) 1238 i=i+1 1239c ipos1=(int_mb(k_offset_alpha+g2b-1)+igi2-1)*nbf+tot_azone1_sh 1240c & +inu1 1241 ipos1=(int_mb(k_offset_alpha+g2b-1)+igi2-1)*nbf+tot_azone2_sh 1242 & +imu1 1243 dbl_mb(k_2g2z+i-1)=dbl_mb(k_movecs_orb+ipos1-1) 1244 enddo 1245 enddo 1246c v(g4b,g3b,inu,imu) C(imu,g2b) 1247 if (.not.ma_push_get(mt_dbl,size_1a3m,'1a3m', 1248 1 l_1a3m_sort,k_1a3m_sort)) 1249 1 call errquit('tce_4ind: step1_1',0,MA_ERR) 1250c dgemm 1251 call dgemm('N','N', 1252 & int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1253 & *nalength(azone1), 1254 & int_mb(k_range_alpha+g2b-1), 1255 & nalength(azone2), 1256 & 1.0d0,dbl_mb(k_2a2m), 1257 & int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1258 & *nalength(azone1), 1259 & dbl_mb(k_2g2z),nalength(azone2),0.0d0, 1260 & dbl_mb(k_1a3m_sort), 1261 & int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1262 & *nalength(azone1)) 1263c transpose 1264 i=0 1265 do inu1=1,nalength(azone1) 1266 do igi2=1,int_mb(k_range_alpha+g2b-1) 1267 do igi3=1,int_mb(k_range_alpha+g3b-1) 1268 do igi4=1,int_mb(k_range_alpha+g4b-1) 1269 i=i+1 1270 ipos1=igi4+int_mb(k_range_alpha+g4b-1)*( 1271 & igi3-1+int_mb(k_range_alpha+g3b-1)*( 1272 & inu1-1+nalength(azone1)*(igi2-1))) 1273 dbl_mb(k_1a3m+i-1)=dbl_mb(k_1a3m+i-1)+ 1274 & dbl_mb(k_1a3m_sort+ipos1-1) 1275 enddo 1276 enddo 1277 enddo 1278 enddo 1279c 1280 if (.not.ma_pop_stack(l_1a3m_sort)) 1281 1 call errquit('tce_mo2e_trans_zones: abc-MA problem',15,MA_ERR) 1282 if (.not.ma_pop_stack(l_2g2z)) 1283 1 call errquit('tce_mo2e_trans_zones: uu l_2g2z',15,MA_ERR) 1284 if (.not.ma_pop_stack(l_2a2m)) 1285 1 call errquit('tce_mo2e_trans_zones: uu l_2g2z',15,MA_ERR) 1286 ENDDO ! PART2 ends up here ------------------------------------------ 1287c 1288c 1289c 1290c write to file 1291 key_1a3m= azone1-1+atpart*(g2b-1+(noa+nva)* 1292 & (g4b-1+(noa+nva)*(g3b-1))) 1293 call tce_hash(int_mb(k_1a3m_offset),key_1a3m,offset_1a3m) 1294 if(idiskl) then !--------- 1295 if (sf_write(i_to_x,dble(bytes)*dble(offset_1a3m), 1296 r dble(bytes)*dble(size_1a3m),dbl_mb(k_1a3m),request).ne.0) 1297 2 call errquit('zones put: sf problem21-b',1,DISK_ERR) 1298 if (sf_wait(request).ne.0) 1299 1 call errquit('zones put: sf problem31-b',2,DISK_ERR) 1300 else 1301 call ga_put(i_to_x,offset_1a3m+1,offset_1a3m+size_1a3m,1,1, 1302 1 dbl_mb(k_1a3m),1) 1303 end if ! ------------- 1304c 1305 if (.not.ma_pop_stack(l_1a3m)) 1306 1 call errquit('tce_mo2e_trans_zones: abc-MA problem',15,MA_ERR) 1307c 1308 next = NXTASK(nprocs, 1) 1309 END IF 1310 count = count + 1 1311 ENDDO 1312 ENDDO 1313 ENDDO 1314 ENDDO 1315c 1A3P fully done here 1316 call ga_sync() 1317 next = nxtask(-nprocs,1) 1318c 1319c 1320c 1321c 1322c 1323 if(nodezero) then 1324 write(6,*)'STEP5 4index' 1325 call util_flush(6) 1326 end if 1327c 1328c 1329c 1330c 1331c 4M starts here 1332c 1333c 1334ccx if(idiskl) then !---------------- 1335ccx i_from_x=i_to_x 1336ccx else !------------------ 1337ccx i_from_x=i_to_x 1338ccx if(i_from_x.eq.d_max1) then 1339ccx i_to_x=d_max2 1340ccx d_v2=i_to_x 1341c size check-out 1342ccx if(size_max2.lt.size_2e) 1343ccx 1 call errquit('tce_4s: 5a',15,MA_ERR) 1344ccx else 1345ccx i_to_x=d_max1 1346ccx d_v2=i_to_x 1347c size check-out 1348ccx if(size_max1.lt.size_2e) 1349ccx 1 call errquit('tce_4s: 5b',15,MA_ERR) 1350ccx end if 1351ccx end if !---------------------- 1352c 1353c DISK CHANGES =============== 4M 1354c destroying handle_2a2m 1355c closing/reopening handle_1a3m file 1356c opening handle_4m file 1357 if(idiskl) then 1358c 1359 if (sf_destroy(handle_2a2m).ne.0) 1360 1 call errquit('tce_sf_destroy: sf problem handle_2a2m',15,MA_ERR) 1361c 1362 if(sf_rwtor(handle_1a3m).ne.0) 1363 2 call errquit('4-index: sf_rwtor 1a3m',0,DISK_ERR) 1364c 1365 i_from_x=handle_1a3m 1366 end if 1367c ============================ 1368c 1369c 1370 nprocs = GA_NNODES() 1371 count = 0 1372 next = NXTASK(nprocs, 1) 1373c 1374 DO g1b = 1,noa+nva !nu 1375 DO g2b = g1b,noa+nva ! 1376 DO g3b = 1,noa+nva !k 1377 DO g4b = g3b,noa+nva !l 1378 IF (next.eq.count) THEN 1379 IF (int_mb(k_spin_alpha+g3b-1)+int_mb(k_spin_alpha+g4b-1).eq. 1380 &int_mb(k_spin_alpha+g1b-1)+int_mb(k_spin_alpha+g2b-1)) THEN 1381 IF (ieor(int_mb(k_sym_alpha+g3b-1),ieor(int_mb(k_sym_alpha+g4b-1), 1382 & ieor(int_mb(k_sym_alpha+g1b-1),int_mb(k_sym_alpha+g2b-1)))) 1383 & .eq. irrep_v) THEN 1384c reversed order 1385 ICOL=INDEX_PAIR(g4b,g3b) 1386 IROW=INDEX_PAIR(g2b,g1b) 1387 IF(IROW.GE.ICOL) THEN 1388 IRES=INDEX_PAIR(IROW,ICOL) 1389 size_4m=int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1390 1 *int_mb(k_range_alpha+g2b-1)*int_mb(k_range_alpha+g1b-1) 1391 if (.not.ma_push_get(mt_dbl,size_4m,'4m',l_4m,k_4m)) 1392 1 call errquit('tce_4ind: step1_1',0,MA_ERR) 1393 call dfill(size_4m, 0.0d0, dbl_mb(k_4m), 1) 1394c 1395 call a1_offset(g1b,g2b,g4b,size_ic, 1396 & l_a1_offset,k_a1_offset) 1397 length_a1=int_mb(k_a1_offset) 1398 DO ii=1,length_a1 ! ii---loop 1399 a1_ini = int_mb(k_a1_offset+ii) 1400 a1_fin = int_mb(k_a1_offset+length_a1+ii) 1401 key_ini = a1_ini-1+atpart*(g4b-1+(noa+nva)* 1402 & (g2b-1+(noa+nva)*(g1b-1))) 1403 call tce_hash(int_mb(k_1a3m_offset),key_ini,offset_bl1a3m) 1404c 1405 size_bl1a3m = int_mb(k_a1_offset+2*length_a1+ii) 1406 if(idiskl) then ! ------- 1407 if (sf_read(i_from_x,dble(bytes)*dble(offset_bl1a3m), 1408 2 dble(bytes)*dble(size_bl1a3m), 1409 2 dbl_mb(k_integral),request).ne.0) 1410 2 call errquit('zones put: sf problem21-b',1,DISK_ERR) 1411 if (sf_wait(request).ne.0) 1412 1 call errquit('zones put: sf problem31-b',2,DISK_ERR) 1413 else 1414 call ga_get(i_from_x,offset_bl1a3m+1,offset_bl1a3m+size_bl1a3m, 1415 1 1,1,dbl_mb(k_integral),1) 1416 end if ! ------- 1417 offset_ismall=0 1418 do i=a1_ini,a1_fin ! over small i ------------------ 1419c C(azone1 g3b) 1420 j=0 1421 tot_azone1_sh=tot_zone(i) 1422 do igi3=1,int_mb(k_range_alpha+g3b-1) 1423 do inu1 = 1,nalength(i) 1424 j=j+1 1425 ipos1=(int_mb(k_offset_alpha+g3b-1)+igi3-1)*nbf+tot_azone1_sh 1426 & +inu1 1427 dbl_mb(k_coeff+j-1)=dbl_mb(k_movecs_orb+ipos1-1) 1428 enddo 1429 enddo 1430c v(g2b,g1b,g4b,inu)*C(inu g3b) 1431 call dgemm('N','N', 1432 & int_mb(k_range_alpha+g2b-1)*int_mb(k_range_alpha+g1b-1) 1433 & *int_mb(k_range_alpha+g4b-1), 1434 & int_mb(k_range_alpha+g3b-1), 1435 & nalength(i), 1436 & 1.0d0,dbl_mb(k_integral+offset_ismall), 1437 & int_mb(k_range_alpha+g2b-1)*int_mb(k_range_alpha+g1b-1) 1438 & *int_mb(k_range_alpha+g4b-1), 1439 & dbl_mb(k_coeff),nalength(i),1.0d0, 1440 & dbl_mb(k_4m), 1441 & int_mb(k_range_alpha+g2b-1)*int_mb(k_range_alpha+g1b-1) 1442 & *int_mb(k_range_alpha+g4b-1)) 1443c 1444 offset_ismall=offset_ismall+nalength(i) 1445 1 *int_mb(k_range_alpha+g4b-1) 1446 1 *int_mb(k_range_alpha+g1b-1)*int_mb(k_range_alpha+g2b-1) 1447c 1448 enddo ! over small i ---------------------- 1449 ENDDO ! ii---loop 1450 if (.not.ma_pop_stack(l_a1_offset)) 1451 1 call errquit('tce_mo2e_trans_zones: uu l_2g2z',15,MA_ERR) 1452c 1453c zapis 1454c 1455 key_4m=ires 1456 call tce_hash_n(int_mb(kax_v2_alpha_offset),key_4m,offset_4m) 1457 call put_block(d_v2,dbl_mb(k_4m),size_4m, 1458 & offset_4m) 1459c 1460 if (.not.ma_pop_stack(l_4m)) 1461 1 call errquit('tce_mo2e_trans_zones: abc-MA problem',15,MA_ERR) 1462 1463 END IF 1464 END IF 1465 END IF 1466 next = NXTASK(nprocs, 1) 1467 END IF 1468 count = count + 1 1469 ENDDO 1470 ENDDO 1471 ENDDO 1472 ENDDO 1473c 1474 call ga_sync() 1475 next = nxtask(-nprocs,1) 1476c 1477c 1478c 1479c 1480c 1481c 1482c 1483c 1484c 1485 if(nodezero) then 1486 write(6,*)'DONE --- DONE ---- DONE ---- DONE' 1487 end if 1488c 1489 if (.not.ma_pop_stack(l_work2)) 1490 1 call errquit('tce_ao2e: MA problem',14,MA_ERR) 1491 if (.not.ma_pop_stack(l_work1)) 1492 1 call errquit('tce_ao2e: MA problem',15,MA_ERR) 1493c 1494 if (.not.ma_pop_stack(l_movecs_orb)) 1495 1 call errquit('tce_ao2e: MA problem',15,MA_ERR) 1496c 1497c 1498c DISK CHANGES =============== AFTER 4M 1499c destroying handle_1a3m 1500 if(idiskl) then 1501c 1502 if (sf_destroy(handle_1a3m).ne.0) 1503 1 call errquit('tce_sf_destroy: sf problem handle_1a3m',15,MA_ERR) 1504c 1505 end if 1506c ============================ 1507c 1508c 1509 if(idiskl) then 1510c if (.not.sf_destroy(d_max2)) 1511c 1 call errquit('tce_sf_destroy22: sf problem',15,MA_ERR) 1512c if (parallel) call ga_sync() 1513c if (.not.sf_destroy(d_max1)) 1514c 1 call errquit('tce_sf_destroy33: sf problem',15,MA_ERR) 1515c if (parallel) call ga_sync() 1516 else 1517 call deletefile(i_from_x) 1518 call ga_sync() 1519cccx call deletefile(d_4af) 1520cccx call ga_sync() 1521 end if 1522c 1523 if (.not.ma_pop_stack(l_1a3m_offset)) 1524 1 call errquit('tce_off_1a3m: MA problem',15,MA_ERR) 1525c 1526 if (.not.ma_pop_stack(l_2a2m_offset)) 1527 1 call errquit('tce_off_2a2m: MA problem',15,MA_ERR) 1528c 1529 if (.not.ma_pop_stack(l_3a1m_offset)) 1530 1 call errquit('tce_off_3a1m: MA problem',15,MA_ERR) 1531c 1532 if (.not.ma_pop_stack(l_4af_offset)) 1533 1 call errquit('tce_off_4a: MA problem',15,MA_ERR) 1534c 1535c 1536c closing two large local files 1537c 1538 if (.not.ma_pop_stack(l_coeff)) 1539 1 call errquit('tce_off_4a: MA problem',15,MA_ERR) 1540c 1541 if (.not.ma_pop_stack(l_integral)) 1542 1 call errquit('tce_off_4a: MA problem',15,MA_ERR) 1543c 1544 1545c 1546 call ga_sync() 1547c *** debug *** 1548c 800 format('DGEMM1 MAX',i5,2x,3f15.5) 1549c 801 format('DGEMM2 ',i5,2x,3f15.5) 1550 8153 format('SF_READ ERROR CODE = ',2x,i10) 1551 8154 format('SF_READ_WAIT ERROR CODE = ',2x,i10) 1552 8156 format('SF_WRITE ERROR CODE = ',2x,i10) 1553 8157 format('SF_WRITE_WAIT ERROR CODE = ',2x,i10) 1554 8155 format('FU',i6,2x,3i20) 1555 9000 format('PART1',i4,1x,'Cpu wall ',2(f17.12,1x),3x,'g4b g3b',2i5) 1556c 9001 format('PART2',i4,1x,'Cpu wall ',2(f17.12,1x),3x,'g4b g3b',2i5) 1557c 9003 format('PART1-4a',i4,1x,'Cpu wall ',2(f17.12,1x)) 1558c 9004 format('PART1-2g2z',i4,1x,'Cpu wall ',2(f17.12,1x)) 1559c 9005 format('PART1-dgemm',i4,1x,'Cpu wall ',2(f17.12,1x)) 1560c 9010 format(' P1-mnrs',i3,1x,2i5,1x,2i5,1x,'Cpu wall ',2(f17.12,1x)) 1561 555 format('atom loop ',2x,i5,3x,2i5,3x,2i5,i12) 1562 556 format('atom time',2x,i5,3x,2i5,3x,2i5,'Cpu wall ',2(f12.7,1x)) 1563 777 format('main do loop ',2x,i5,3x,2i5,3x,2i5) 1564 775 format('main loop step1 ',2x,i5,3x,2i5,3x,2i5) 1565 776 format('main loop step2 ',2x,i5,3x,2i5,3x,2i5) 1566 778 format('PART1',2x,i5,3x,2i5,3x,2i5,2x,'Cpu wall ',2(f17.12,1x)) 1567 779 format('PART2',2x,i5,3x,2i5,3x,2i5,2x,'Cpu wall ',2(f17.12,1x)) 1568 780 format('ADD BLOCK',2x,i5,3x,2i5,3x,2i5,2x,'Cpu wall ', 1569 & 2(f17.12,1x)) 1570c ************* 1571c 1572 RETURN 1573 END 1574c 1575c 1576c 1577 SUBROUTINE tce_4a_offsetx2(l_a_offset,k_a_offset,size) 1578C $Id$ 1579C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 1580C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 1581 IMPLICIT NONE 1582#include "global.fh" 1583#include "mafdecls.fh" 1584#include "sym.fh" 1585#include "errquit.fh" 1586#include "tce.fh" 1587#include "tce_main.fh" 1588 INTEGER l_a_offset 1589 INTEGER k_a_offset 1590 INTEGER size 1591 INTEGER length 1592 INTEGER addr 1593 INTEGER mu 1594 INTEGER nu 1595 INTEGER rho 1596 INTEGER sigma 1597 INTEGER azone1,azone2,azone3,azone4 1598 length = 0 1599 DO azone1 = 1,atpart !nu 1600 DO azone2 = azone1,atpart !mu 1601 DO azone3 = 1,atpart !sigma 1602 DO azone4 = azone3,atpart !rho 1603 length = length + 1 1604 END DO 1605 END DO 1606 END DO 1607 END DO 1608 IF (.not.MA_PUSH_GET(mt_int,2*length+1,'noname',l_a_offset,k_a_off 1609 &set)) CALL ERRQUIT('tce_t2_offset',0,MA_ERR) 1610 int_mb(k_a_offset) = length 1611 addr = 0 1612 size = 0 1613 DO azone1 = 1,atpart !nu 1614 DO azone2 = azone1,atpart !mu 1615 DO azone3 = 1,atpart !sigma 1616 DO azone4 = azone3,atpart !rho 1617 addr = addr + 1 1618 int_mb(k_a_offset+addr) = azone4 - 1 + atpart * (azone3 - 1 + 1619 & atpart * (azone2 - 1 + atpart * (azone1 - 1))) 1620 int_mb(k_a_offset+length+addr) = size 1621 size = size + nalength(azone1) * nalength(azone2) * 1622 & nalength(azone3) * nalength(azone4) 1623 END DO 1624 END DO 1625 END DO 1626 END DO 1627 RETURN 1628 END 1629c 1630c 1631 SUBROUTINE tce_3a1m_offsetx2(l_a_offset,k_a_offset,size) 1632C $Id$ 1633C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 1634C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 1635 IMPLICIT NONE 1636#include "global.fh" 1637#include "mafdecls.fh" 1638#include "sym.fh" 1639#include "errquit.fh" 1640#include "tce.fh" 1641#include "tce_main.fh" 1642 INTEGER l_a_offset 1643 INTEGER k_a_offset 1644 INTEGER size 1645 INTEGER length 1646 INTEGER addr 1647 INTEGER mu 1648 INTEGER nu 1649 INTEGER rho 1650 INTEGER sigma 1651 INTEGER azone1,azone2,azone3 1652 INTEGER g2b 1653 length = 0 1654 DO azone1 = 1,atpart !nu 1655 DO azone2 = azone1,atpart !mu 1656 DO g2b = 1,noa+nva !g2b 1657 DO azone3 = 1,atpart !sigma 1658 length = length + 1 1659 END DO 1660 END DO 1661 END DO 1662 END DO 1663 IF (.not.MA_PUSH_GET(mt_int,2*length+1,'noname',l_a_offset,k_a_off 1664 &set)) CALL ERRQUIT('tce_t2_offset',0,MA_ERR) 1665 int_mb(k_a_offset) = length 1666 addr = 0 1667 size = 0 1668 DO azone1 = 1,atpart !nu 1669 DO azone2 = azone1,atpart !mu 1670 DO g2b = 1,noa+nva !g2b 1671 DO azone3 = 1,atpart !sigma 1672 addr = addr + 1 1673 int_mb(k_a_offset+addr) = azone3-1+atpart*(g2b-1+(noa+nva)* 1674 & (azone2 - 1 + atpart * (azone1 - 1))) 1675 int_mb(k_a_offset+length+addr) = size 1676 size = size + nalength(azone1) * nalength(azone2) * 1677 & nalength(azone3) * int_mb(k_range_alpha+g2b-1) 1678 END DO 1679 END DO 1680 END DO 1681 END DO 1682 RETURN 1683 END 1684c 1685c 1686c 1687c 1688 SUBROUTINE tce_2a2m_offsetx2(l_a_offset,k_a_offset,size) 1689C $Id$ 1690C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 1691C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 1692 IMPLICIT NONE 1693#include "global.fh" 1694#include "mafdecls.fh" 1695#include "sym.fh" 1696#include "errquit.fh" 1697#include "tce.fh" 1698#include "tce_main.fh" 1699 INTEGER l_a_offset 1700 INTEGER k_a_offset 1701 INTEGER size 1702 INTEGER length 1703 INTEGER addr 1704 INTEGER g1b,g2b 1705 INTEGER azone1,azone2 1706 length = 0 1707 DO g1b = 1,noa+nva !k 1708 DO g2b = g1b,noa+nva !l 1709 DO azone1 = 1,atpart !nu 1710 DO azone2 = azone1,atpart !mu 1711 length = length + 1 1712 END DO 1713 END DO 1714 END DO 1715 END DO 1716 IF (.not.MA_PUSH_GET(mt_int,2*length+1,'noname',l_a_offset,k_a_off 1717 &set)) CALL ERRQUIT('tce_t2_offset',0,MA_ERR) 1718 int_mb(k_a_offset) = length 1719 addr = 0 1720 size = 0 1721 DO g1b = 1,noa+nva !k 1722 DO g2b = g1b,noa+nva !l 1723 DO azone1 = 1,atpart !nu 1724 DO azone2 = azone1,atpart !mu 1725 addr = addr + 1 1726 int_mb(k_a_offset+addr) = azone2 - 1 + atpart*( azone1 - 1 + 1727 & atpart*(g2b - 1 + (noa+nva) * (g1b-1))) 1728 int_mb(k_a_offset+length+addr) = size 1729 size = size + 1730 & int_mb(k_range_alpha+g1b-1) * int_mb(k_range_alpha+g2b-1) 1731 & * nalength(azone1) * nalength(azone2) 1732 END DO 1733 END DO 1734 END DO 1735 END DO 1736 RETURN 1737 END 1738c 1739c 1740c 1741 SUBROUTINE tce_1a3m_offsetx2(l_a_offset,k_a_offset,size) 1742C $Id$ 1743C This is a Fortran77 program generated by Tensor Contraction Engine v.1.0 1744C Copyright (c) Battelle & Pacific Northwest National Laboratory (2002) 1745 IMPLICIT NONE 1746#include "global.fh" 1747#include "mafdecls.fh" 1748#include "sym.fh" 1749#include "errquit.fh" 1750#include "tce.fh" 1751#include "tce_main.fh" 1752 INTEGER l_a_offset 1753 INTEGER k_a_offset 1754 INTEGER size 1755 INTEGER length 1756 INTEGER addr 1757 INTEGER g1b,g2b,g4b 1758 INTEGER azone1,azone2 1759 length = 0 1760 DO g1b = 1,noa+nva !k 1761 DO g2b = g1b,noa+nva !l 1762 DO g4b = 1,noa+nva !i 1763 DO azone1 = 1,atpart !nu 1764 length = length + 1 1765 END DO 1766 END DO 1767 END DO 1768 END DO 1769 IF (.not.MA_PUSH_GET(mt_int,2*length+1,'noname',l_a_offset,k_a_off 1770 &set)) CALL ERRQUIT('tce_t2_offset',0,MA_ERR) 1771 int_mb(k_a_offset) = length 1772 addr = 0 1773 size = 0 1774 DO g1b = 1,noa+nva !k 1775 DO g2b = g1b,noa+nva !l 1776 DO g4b = 1,noa+nva !i 1777 DO azone1 = 1,atpart !nu 1778 addr = addr + 1 1779 int_mb(k_a_offset+addr) = azone1-1+atpart*(g4b-1+(noa+nva)* 1780 & (g2b-1+(noa+nva)*(g1b-1))) 1781 int_mb(k_a_offset+length+addr) = size 1782 size = size + 1783 & int_mb(k_range_alpha+g1b-1) * int_mb(k_range_alpha+g2b-1) 1784 & * int_mb(k_range_alpha+g4b-1) * nalength(azone1) 1785 END DO 1786 END DO 1787 END DO 1788 END DO 1789 RETURN 1790 END 1791c 1792 subroutine a4_offset(azone1,azone2,azone3,size_ic, 1793 & l_a4_offset,k_a4_offset) 1794 IMPLICIT NONE 1795#include "global.fh" 1796#include "mafdecls.fh" 1797#include "sym.fh" 1798#include "errquit.fh" 1799#include "tce.fh" 1800#include "tce_main.fh" 1801 INTEGER azone1,azone2,azone3,size_ic,l_a4_offset,k_a4_offset 1802 INTEGER i,j,k,size_int,length_a4,azone4,a4_ini,a4_fin,offset_a4 1803 length_a4=0 1804 size_int=0 1805 do azone4=azone3,atpart 1806 size_int=size_int+nalength(azone1)*nalength(azone2)* 1807 1 nalength(azone3)*nalength(azone4) 1808 if(size_int.gt.size_ic) then 1809 length_a4=length_a4+1 1810 size_int= nalength(azone1)*nalength(azone2)* 1811 1 nalength(azone3)*nalength(azone4) 1812 end if 1813 enddo 1814 if(size_int.ne.0) length_a4=length_a4+1 1815 IF (.not.MA_PUSH_GET(mt_int,3*length_a4+1,'noname', 1816 &l_a4_offset,k_a4_offset)) 1817 & CALL ERRQUIT('tce_a4_offset',0,MA_ERR) 1818 int_mb(k_a4_offset) = length_a4 1819c 1820 i=0 1821 a4_ini=azone3 1822 size_int=0 1823 do azone4=azone3,atpart 1824 size_int=size_int+nalength(azone1)*nalength(azone2)* 1825 1 nalength(azone3)*nalength(azone4) 1826 if(size_int.gt.size_ic) then 1827 i=i+1 1828 int_mb(k_a4_offset+i)=a4_ini 1829 int_mb(k_a4_offset+length_a4+i)=azone4-1 1830 int_mb(k_a4_offset+2*length_a4+i)=size_int- 1831 1 nalength(azone1)*nalength(azone2)* 1832 1 nalength(azone3)*nalength(azone4) 1833 size_int= nalength(azone1)*nalength(azone2)* 1834 1 nalength(azone3)*nalength(azone4) 1835 a4_ini=azone4 1836 end if 1837 enddo 1838 if(size_int.ne.0) then 1839 i=i+1 1840 int_mb(k_a4_offset+i)=a4_ini 1841 int_mb(k_a4_offset+length_a4+i)=atpart 1842 int_mb(k_a4_offset+2*length_a4+i)=size_int 1843 end if 1844 if(i.ne.length_a4) 1845 & CALL ERRQUIT('tce_a4_off problems',0,MA_ERR) 1846 return 1847 end 1848c 1849c 1850c 1851 subroutine a3_offset(azone1,azone2,g3b,size_ic, 1852 & l_a3_offset,k_a3_offset) 1853 IMPLICIT NONE 1854#include "global.fh" 1855#include "mafdecls.fh" 1856#include "sym.fh" 1857#include "errquit.fh" 1858#include "tce.fh" 1859#include "tce_main.fh" 1860 INTEGER azone1,azone2,g3b,size_ic,l_a3_offset,k_a3_offset 1861 INTEGER i,j,k,size_int,length_a3,azone3,a3_ini,a3_fin,offset_a3 1862 length_a3=0 1863 size_int=0 1864 do azone3=1,atpart 1865 size_int=size_int+nalength(azone1)*nalength(azone2)* 1866 1 nalength(azone3)*int_mb(k_range_alpha+g3b-1) 1867 if(size_int.gt.size_ic) then 1868 length_a3=length_a3+1 1869 size_int= nalength(azone1)*nalength(azone2)* 1870 1 nalength(azone3)*int_mb(k_range_alpha+g3b-1) 1871 end if 1872 enddo 1873 if(size_int.ne.0) length_a3=length_a3+1 1874 IF (.not.MA_PUSH_GET(mt_int,3*length_a3+1,'noname', 1875 &l_a3_offset,k_a3_offset)) 1876 & CALL ERRQUIT('tce_a3_offset',0,MA_ERR) 1877 int_mb(k_a3_offset) = length_a3 1878c 1879 i=0 1880 a3_ini=1 1881 size_int=0 1882 do azone3=1,atpart 1883 size_int=size_int+nalength(azone1)*nalength(azone2)* 1884 1 nalength(azone3)*int_mb(k_range_alpha+g3b-1) 1885 if(size_int.gt.size_ic) then 1886 i=i+1 1887 int_mb(k_a3_offset+i)=a3_ini 1888 int_mb(k_a3_offset+length_a3+i)=azone3-1 1889 int_mb(k_a3_offset+2*length_a3+i)=size_int- 1890 1 nalength(azone1)*nalength(azone2)* 1891 1 nalength(azone3)*int_mb(k_range_alpha+g3b-1) 1892 size_int= nalength(azone1)*nalength(azone2)* 1893 1 nalength(azone3)*int_mb(k_range_alpha+g3b-1) 1894 a3_ini=azone3 1895 end if 1896 enddo 1897 if(size_int.ne.0) then 1898 i=i+1 1899 int_mb(k_a3_offset+i)=a3_ini 1900 int_mb(k_a3_offset+length_a3+i)=atpart 1901 int_mb(k_a3_offset+2*length_a3+i)=size_int 1902 end if 1903 if(i.ne.length_a3) 1904 & CALL ERRQUIT('tce_a3_off problems',0,MA_ERR) 1905 return 1906 end 1907c 1908c 1909c 1910 subroutine a2_offset(azone1,g4b,g3b,size_ic, 1911 & l_a2_offset,k_a2_offset) 1912 IMPLICIT NONE 1913#include "global.fh" 1914#include "mafdecls.fh" 1915#include "sym.fh" 1916#include "errquit.fh" 1917#include "tce.fh" 1918#include "tce_main.fh" 1919 INTEGER azone1,g4b,g3b,size_ic,l_a2_offset,k_a2_offset 1920 INTEGER i,j,k,size_int,length_a2,azone2,a2_ini,a2_fin,offset_a2 1921 length_a2=0 1922 size_int=0 1923 do azone2=azone1,atpart 1924 size_int=size_int+nalength(azone1)*nalength(azone2)* 1925 1 int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1926 if(size_int.gt.size_ic) then 1927 length_a2=length_a2+1 1928 size_int= nalength(azone1)*nalength(azone2)* 1929 1 int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1930 end if 1931 enddo 1932 if(size_int.ne.0) length_a2=length_a2+1 1933 IF (.not.MA_PUSH_GET(mt_int,3*length_a2+1,'noname', 1934 &l_a2_offset,k_a2_offset)) 1935 & CALL ERRQUIT('tce_a2_offset',0,MA_ERR) 1936 int_mb(k_a2_offset) = length_a2 1937c 1938 i=0 1939 a2_ini=azone1 1940 size_int=0 1941 do azone2=azone1,atpart 1942 size_int=size_int+nalength(azone1)*nalength(azone2)* 1943 1 int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1944 if(size_int.gt.size_ic) then 1945 i=i+1 1946 int_mb(k_a2_offset+i)=a2_ini 1947 int_mb(k_a2_offset+length_a2+i)=azone2-1 1948 int_mb(k_a2_offset+2*length_a2+i)=size_int- 1949 1 nalength(azone1)*nalength(azone2)* 1950 1 int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1951 size_int= nalength(azone1)*nalength(azone2)* 1952 1 int_mb(k_range_alpha+g4b-1)*int_mb(k_range_alpha+g3b-1) 1953 a2_ini=azone2 1954 end if 1955 enddo 1956 if(size_int.ne.0) then 1957 i=i+1 1958 int_mb(k_a2_offset+i)=a2_ini 1959 int_mb(k_a2_offset+length_a2+i)=atpart 1960 int_mb(k_a2_offset+2*length_a2+i)=size_int 1961 end if 1962 if(i.ne.length_a2) 1963 & CALL ERRQUIT('tce_a2_off problems',0,MA_ERR) 1964 return 1965 end 1966c 1967c 1968c 1969 subroutine a1_offset(g1b,g2b,g3b,size_ic, 1970 & l_a1_offset,k_a1_offset) 1971 IMPLICIT NONE 1972#include "global.fh" 1973#include "mafdecls.fh" 1974#include "sym.fh" 1975#include "errquit.fh" 1976#include "tce.fh" 1977#include "tce_main.fh" 1978 INTEGER g1b,g2b,g3b,size_ic,l_a1_offset,k_a1_offset 1979 INTEGER i,j,k,size_int,length_a1,azone1,a1_ini,a1_fin,offset_a1 1980 length_a1=0 1981 size_int=0 1982 do azone1=1,atpart 1983 size_int=size_int+nalength(azone1)*int_mb(k_range_alpha+g1b-1)* 1984 1 int_mb(k_range_alpha+g2b-1)*int_mb(k_range_alpha+g3b-1) 1985 if(size_int.gt.size_ic) then 1986 length_a1=length_a1+1 1987 size_int= nalength(azone1)*int_mb(k_range_alpha+g1b-1)* 1988 1 int_mb(k_range_alpha+g2b-1)*int_mb(k_range_alpha+g3b-1) 1989 end if 1990 enddo 1991 if(size_int.ne.0) length_a1=length_a1+1 1992 IF (.not.MA_PUSH_GET(mt_int,3*length_a1+1,'noname', 1993 &l_a1_offset,k_a1_offset)) 1994 & CALL ERRQUIT('tce_a1_offset',0,MA_ERR) 1995 int_mb(k_a1_offset) = length_a1 1996c 1997 i=0 1998 a1_ini=1 1999 size_int=0 2000 do azone1=1,atpart 2001 size_int=size_int+nalength(azone1)*int_mb(k_range_alpha+g1b-1)* 2002 1 int_mb(k_range_alpha+g2b-1)*int_mb(k_range_alpha+g3b-1) 2003 if(size_int.gt.size_ic) then 2004 i=i+1 2005 int_mb(k_a1_offset+i)=a1_ini 2006 int_mb(k_a1_offset+length_a1+i)=azone1-1 2007 int_mb(k_a1_offset+2*length_a1+i)=size_int- 2008 1 nalength(azone1)*int_mb(k_range_alpha+g1b-1)* 2009 1 int_mb(k_range_alpha+g2b-1)*int_mb(k_range_alpha+g3b-1) 2010 size_int=nalength(azone1)*int_mb(k_range_alpha+g1b-1)* 2011 1 int_mb(k_range_alpha+g2b-1)*int_mb(k_range_alpha+g3b-1) 2012 a1_ini=azone1 2013 end if 2014 enddo 2015 if(size_int.ne.0) then 2016 i=i+1 2017 int_mb(k_a1_offset+i)=a1_ini 2018 int_mb(k_a1_offset+length_a1+i)=atpart 2019 int_mb(k_a1_offset+2*length_a1+i)=size_int 2020 end if 2021 if(i.ne.length_a1) 2022 & CALL ERRQUIT('tce_a1_off problems',0,MA_ERR) 2023 return 2024 end 2025