1 Subroutine hf1(Axyz,Aprims,Acoefs,NPA,NCA,La, 2 & Bxyz,Bprims,Bcoefs,NPB,NCB,Lb, 3 & Cxyz,zan,exinv,ncenters, 4 & bO2I,bKEI,bNAI,Nints,O2I,KEI,NAI,canAB, 5 & DryRun,W0,maxW0) 6c $Id$ 7 Implicit none 8 integer NPA,NCA,La,NPB,NCB,Lb,ncenters,Nints,maxW0 9 logical O2I,KEI,NAI,canAB 10 logical GenCon,DryRun 11 12c--> Cartesian Coordinates, Primitives & Contraction Coefficients 13 14 double precision Axyz(3),Aprims(NPA),Acoefs(NPA,NCA) 15 double precision Bxyz(3),Bprims(NPB),Bcoefs(NPB,NCB) 16 17c--> Nuclear Cartesian Coordinates, Charges and Inverse Exponents 18 19 double precision Cxyz(3,ncenters),zan(ncenters),exinv(ncenters) 20 21c--> Blocks of Overlap, Kinetic Energy & Nuclear Attraction Integrals 22 23 double precision bO2I(Nints),bKEI(Nints),bNAI(Nints) 24 25c--> Scratch Space. 26 27 double precision W0(maxW0) 28 29c--> Local variables 30 31 integer MXD,NCP,NPP,La2,Lb2,Li,Lp,Lp3,MaxMem,nintlocal,lprod,nd 32 integer i_ALPHAp,i_IPAIRp,i_ESp,i_left,i_right,i_top,i_Ep,i_ptr, 33 & i_pf,i_prim_ints,i_half_ints,i_Ti,i_R0,i_IJK,i_P,i_RS,i_PC, 34 & i_ff,i_Rj 35#if defined(INTDEBUG) 36 integer jjjj,iiii 37#endif 38#include "mafdecls.fh" 39#include "stdio.fh" 40#include "errquit.fh" 41c 42c Compute the overlap, kinetic energy, and nuclear attraction integrals for 43c two shells of contracted Gaussians functions. This driver is NOT capable of 44c evaluating integral derivatives. 45c 46c****************************************************************************** 47#if defined(INTDEBUG) 48 if (.not.dryrun) then 49 write(LuOut,*)' inside hf1 ' 50 write(LuOut,*)' npa,nca,la = ',npa,nca,la 51 write(LuOut,*)' npb,ncb,lb = ',npb,ncb,lb 52 write(LuOut,*)' ncenters = ',ncenters 53 write(LuOut,*)' NINTS = ',nints 54 write(LuOut,*)' maxW0 = ',maxw0 55 write(LuOut,*)' <canAB:DryRun>-<',canab,':',dryrun,'>' 56 write(LuOut,*)' <o2i:kei:nai>-<',o2i,':',kei,':',nai,'>' 57 write(6,'(a,3(2x,1pd20.10))')' Axyz =',Axyz 58 write(6,'(a,3(2x,1pd20.10))')' Bxyz =',Bxyz 59 write(6,'(a,100(3(2x,1pd20.10/)))')' Cxyz =',Cxyz 60 do jjjj = 1,nca 61 do iiii = 1,npa 62 write(6,'(a,i3,a,2(2x,1pd20.10))') 63 & 'Aprims:Acoeffs:(',iiii,') =',Aprims(iiii), 64 & Acoefs(iiii,jjjj) 65 enddo 66 enddo 67 do jjjj = 1,ncb 68 do iiii = 1,npb 69 write(6,'(a,i3,a,2(2x,1pd20.10))') 70 & 'Bprims:Bcoeffs:(',iiii,') =',Bprims(iiii), 71 & Bcoefs(iiii,jjjj) 72 enddo 73 enddo 74 endif 75 if (.not.dryrun) then 76 write(LuOut,*)' Li/Lj',La,'/',Lb 77 write(LuOut,*)' i_ngen ',nca 78 write(LuOut,*)' j_ngen ',ncb 79 write(LuOut,*)' int_hf1: lstv',Nints 80 write(LuOut,*)' int_hf1: lscr',maxW0 81 call dcopy(maxW0,0d0,0,W0,1) 82 if (O2I) call dcopy(Nints,0d0,0,bO2I,1) 83 if (KEI) call dcopy(Nints,0d0,0,bKEI,1) 84 if (NAI) call dcopy(Nints,0d0,0,bNAI,1) 85*debug_ma: call MA_summarize_allocated_blocks() 86*debug_ma: write(LuOut,*)' int_hf1: ma verify 1-b4' 87*debug_ma: status = ma_verify_allocator_stuff() 88*debug_ma: write(LuOut,*)' verstat = ',status 89*debug_ma: write(LuOut,*)' int_hf1: ma verify 1-af' 90 call util_flush(6) 91 endif 92 if (.not.dryrun) then 93 call hf_print_set(1) 94 call hf_print('hf1: a shell',axyz,aprims,acoefs,npa,nca,la) 95 call hf_print('hf1: b shell',bxyz,bprims,bcoefs,npb,ncb,lb) 96 call hf_print_set(0) 97 endif 98#endif 99 100 MXD = 0 101 102c Determine whether general or segmented contraction is used. 103 104 NCP = NCA*NCB 105 GenCon = NCP.ne.1 106 107c To determine all the Hermite expansion coefficients required to evaluate 108c the kinetic energy integrals, increment the angular momenta by one. 109 110 if( KEI )then 111 Li = 1 112 else 113 Li = 0 114 end if 115 116c Define the angular momentum of the overlap distribution. 117 118 Lp = La + Lb 119 120c Increment "Lp" to account for the order of differentiation. 121 122 Lp = Lp + MXD 123 124c Define the accumulated number of angular momentum functions <= Lp. 125 126 Lp3 = ((Lp+1)*(Lp+2)*(Lp+3))/6 127 128c Define the prefactor of the overlap distribution "P". 129 130c Assign pointers to scratch space. 131 132 i_ALPHAp = 1 133 i_IPAIRp = i_ALPHAp + 2*(NPA*NPB) 134 i_left = i_IPAIRp + 2*(NPA*NPB) - 1 135 136 i_ESp = (maxW0+1) - 3*(NPA*NPB) 137 i_right = i_ESp 138 139 MaxMem = 1 ! take care of compiler warnings 140 if( DryRun )then 141 142 MaxMem = i_left + (maxW0 - (i_right-1)) 143 NPP = NPA*NPB 144 145 else if (i_left.ge.i_right) then 146 147 write(LuOut,*) 'HF1: Insufficient scratch space.' 148 write(LuOut,*) ' needed ',i_left+(maxW0-(i_right-1)) 149 write(LuOut,*) ' allocated ',maxW0 150 write(LuOut,*) ' DryRun ',DryRun 151 152 write(LuOut,*) 'From the left ' 153 write(LuOut,*) 'ALPHAp: ',i_ALPHAp 154 write(LuOut,*) 'IPAIRp: ',i_IPAIRp 155 write(LuOut,*) 'From the right ' 156 write(LuOut,*) 'ESp : ',i_ESp 157 158 call errquit('hf1: insufficient memory for hfset',911, MEM_ERR) 159 160 else 161 162 call hfset(Axyz,Aprims,Acoefs,NPA,NCA, 163 & Bxyz,Bprims,Bcoefs,NPB,NCB, 164 & GenCon,W0(i_ALPHAp),W0(i_IPAIRp),W0(i_ESp),NPP) 165 166 end if 167 168 La2 = (La+1)*(La+2)/2 169 Lb2 = (Lb+1)*(Lb+2)/2 170 171c Zero out the integrals. Return if screening test gives no pairs 172 173 if (.not.DryRun) then 174 nintlocal = La2*Lb2*nca*ncb 175 if (O2I) call dcopy(nintlocal,0d0,0,bO2I,1) 176 if (KEI) call dcopy(nintlocal,0d0,0,bKEI,1) 177 if (NAI) call dcopy(nintlocal,0d0,0,bNAI,1) 178 if (NPP.eq.0) return 179 end if 180 181c Define the Hermite linear expansion coefficients. 182 183c Assign pointers to scratch space. 184 185 lprod = ((La+Li)+(Lb+Li)+1)*((La+Li)+1)*((Lb+Li)+1) 186 187C write (LuOut,*) 'before hfmke' 188 189 i_Ep = i_IPAIRp + 2*(NPA*NPB) 190 i_ptr = i_Ep + 3*NPP*(MXD+1)*lprod 191 192 i_prim_ints = i_ptr ! take care of compiler warnings 193 i_half_ints = i_prim_ints + NPP 194 195 if (gencon) then 196 i_prim_ints = i_ptr 197 i_half_ints = i_prim_ints + NPP 198 i_ptr = i_half_ints + NPA*NCB 199 end if 200 201 i_pf = i_ptr 202 i_left = i_pf + 2*NPP - 1 203 204 if( DryRun )then 205 206 MaxMem = max( MaxMem, i_left + (maxW0 - (i_right-1)) ) 207 208 else if( i_left.ge.i_right) then 209 210 write(LuOut,*) 'HF1: Insufficient scratch space.' 211 write(LuOut,*) ' needed ',i_left+(maxW0-(i_right-1)) 212 write(LuOut,*) ' allocated ',maxW0 213 write(LuOut,*) ' DryRun ',DryRun 214 215 write(LuOut,*) 'From the left ' 216 write(LuOut,*) 'ALPHAp: ',i_ALPHAp 217 write(LuOut,*) 'IPAIRp: ',i_IPAIRp 218 write(LuOut,*) 'Ep : ',i_Ep 219 write(LuOut,*) 'pf : ',i_pf 220 write(LuOut,*) 'From the right ' 221 write(LuOut,*) 'ESp : ',i_ESp 222 223 call errquit('hf1: insufficient memory for hfmke',911, MEM_ERR) 224 225 else 226 227 do nd = 0,MXD 228 call hfmke(Axyz,Bxyz,W0(i_ALPHAp),W0(i_ESp),W0(i_Ep),W0(i_pf), 229 & nd,NPP,MXD,La+Li,Lb+Li) 230 end do 231 232 end if 233 234c Compute the 2-center overlap integrals, <a|S|b>. 235 236 if( O2I )then 237 if( .not. DryRun )then 238 if (gencon) then 239 call hf2oi_gc(W0(i_Ep),bO2I,W0(i_prim_ints),W0(i_half_ints), 240 & Acoefs,Bcoefs,W0(i_IPAIRp),NPA,NPB,NCA,NCB,NPP, 241 & La,Lb,La2,Lb2,Li,canAB) 242 else 243 call hf2oi(W0(i_Ep),bO2I,Nints,NPP,La,Lb,Li,canAB) 244 endif 245 end if 246 end if 247 248c Compute kinetic energy integrals, <a|T|b>. 249 250 if (KEI) then 251 252c Assign pointers to scratch space. 253 254 i_Ti = i_ptr 255 i_top = i_Ti + NPP - 1 256 257 if( DryRun )then 258 259 MaxMem = max( MaxMem, i_top ) 260 261 else if( i_top.gt.maxW0 )then 262 263 write(LuOut,*) 'HF1: Insufficient scratch space.' 264 write(LuOut,*) ' needed ',i_top 265 write(LuOut,*) ' allocated ',maxW0 266 write(LuOut,*) ' DryRun ',DryRun 267 268 write(LuOut,*) 'ALPHAp: ',i_ALPHAp 269 write(LuOut,*) 'IPAIRp: ',i_IPAIRp 270 write(LuOut,*) 'Ep : ',i_Ep 271 write(LuOut,*) 'Ti : ',i_Ti 272 273 call errquit('hf1: insufficient memory for hfkei',911, 274 & MEM_ERR) 275 276 else if (gencon) then 277 278 call hfkei_gc(W0(i_ALPHAp),W0(i_Ep),bKEI, 279 & W0(i_prim_ints),W0(i_half_ints),W0(i_Ti), 280 & Acoefs,Bcoefs,W0(i_IPAIRp), 281 & NPA,NPB,NCA,NCB,NPP,La,Lb,La2,Lb2,Li,canAB) 282 else 283 284 call hfkei(W0(i_ALPHAp),W0(i_Ep),bKEI,W0(i_Ti), 285 & Nints,NPP,La,Lb,Li,canAB) 286 287 end if 288 289 end if 290 291c Compute nuclear attraction integrals, <a|V|b>. 292 293 if( NAI )then 294 295c Define the auxiliary function integrals. 296 297c Assign scratch space. 298 299 i_R0 = i_ptr 300 i_IJK = i_R0 + NPP*Lp3 301 i_P = i_IJK + (Lp+1)**3 302 i_RS = i_P + NPP*3 303 i_PC = i_RS + NPP 304 i_ff = i_PC + NPP*3 305 i_Rj = i_ff + NPP*2 306 i_top = i_Rj + NPP*(Lp+1)*Lp3 - 1 307 308 if( DryRun )then 309 310 MaxMem = max( MaxMem, i_top ) 311C write (LuOut,*) MaxMem,maxW0 312 313 else if( i_top.gt.maxW0 .and. .not.DryRun )then 314 315 write(LuOut,*) 'HF1: Insufficient scratch space.' 316 write(LuOut,*) ' needed ',i_top 317 write(LuOut,*) ' allocated ',maxW0 318 write(LuOut,*) ' DryRun ',DryRun 319 320 write(LuOut,*) 'ALPHAp: ',i_ALPHAp 321 write(LuOut,*) 'IPAIRp: ',i_IPAIRp 322 write(LuOut,*) 'Ep : ',i_Ep 323 write(LuOut,*) 'R0 : ',i_R0 324 write(LuOut,*) 'IJK : ',i_IJK 325 write(LuOut,*) 'P : ',i_P 326 write(LuOut,*) 'RS : ',i_RS 327 write(LuOut,*) 'PC : ',i_PC 328 write(LuOut,*) 'ff : ',i_ff 329 write(LuOut,*) 'Rj : ',i_Rj 330 331 call errquit('hf1: insufficient memory for hfmkr',911, 332 & MEM_ERR) 333 334 else 335 336 call hf1mkr(Axyz,Bxyz,Cxyz,zan,exinv,ncenters, 337 & W0(i_ALPHAp),W0(i_P),W0(i_RS),W0(i_PC),W0(i_ff), 338 & W0(i_Rj),W0(i_R0),W0(i_R0),W0(i_IJK), 339 & NPP,Lp,Lp3,.FALSE.) 340 341 if (gencon) then 342 call hfnai_gc(W0(i_Ep),W0(i_R0),W0(i_IJK),bNAI, 343 & W0(i_prim_ints),W0(i_half_ints), 344 & Acoefs,Bcoefs,W0(i_IPAIRp), 345 & NPA,NPB,NCA,NCB,NPP,La,Lb,La2,Lb2,Li,Lp,Lp3,canAB) 346 else 347 call hfnai(W0(i_Ep),W0(i_R0),W0(i_IJK),bNAI, 348 & Nints,NPP,La,Lb,Li,Lp,Lp3,canAB) 349 end if 350 351 end if 352 353 end if 354 355c Return the maximum amount of scratch space required by a "dry run". 356 357 if( DryRun ) maxW0 = MaxMem 358c 359 end 360