1 subroutine hf3pot( 2 & Axyz,Aprims,Acoef,NPA,NCA,La, 3 & Bxyz,Bprims,Bcoef,NPB,NCB,Lb, 4 & Cxyz,Cprims,Ccoef,NPC,NCC,Lc, 5 & b3nai,b3nai_sz,nint, 6 & DryRun,Scr,lscr) 7c 8c $Id$ 9c 10 implicit none 11c:: includes 12#include "stdio.fh" 13#include "mafdecls.fh" 14#include "errquit.fh" 15c:: passed 16 integer La ! [input] angular momentum for center a 17 integer Lb ! [input] angular momentum for center b 18 integer Lc ! [input] angular momentum for center c 19 integer NPA ! [input] number of primitives on center a 20 integer NPB ! [input] number of primitives on center b 21 integer NPC ! [input] number of primitives on center c 22 integer NCA ! [input] number of general contractions on center a 23 integer NCB ! [input] number of general contractions on center b 24 integer NCC ! [input] number of general contractions on center c 25 double precision Axyz(3) ! [input] coordinates for center a 26 double precision Bxyz(3) ! [input] coordinates for center b 27 double precision Cxyz(3) ! [input] coordinates for center c 28 double precision Aprims(NPA) ! [input] primitive exponents on center a 29 double precision Bprims(NPB) ! [input] primitive exponents on center b 30 double precision Cprims(NPC) ! [input] primitive exponents on center c 31 double precision Acoef(NPA,NCA) ! [input] primitive coefficients on center a 32 double precision Bcoef(NPB,NCB) ! [input] primitive coefficients on center b 33 double precision Ccoef(NPC,NCC) ! [input] primitive coefficients on center c 34 integer b3nai_sz ! [input] input size of 3 center nai buffer 35 double precision b3nai(b3nai_sz) ! [output] 3 center nai buffer 36 integer nint ! [output] actual number of integrals computed 37 logical DryRun ! [input] run routine only to test memory size needed 38 integer lscr ! [input] size of scratch array 39 double precision scr(lscr) ! [scratch] scratch array for computations 40c::local 41 logical gencon ! general contraction flag 42 Integer Nabc ! product of primitives 43 Integer maxmem ! memory returned in lscr for dryrun 44 Integer Lg ! total angular momentum 45 Integer Lg3 ! l3 of total angular momentum 46c:: scratch pointers 47 Integer i_top ! top of scratch 48 Integer i_save ! point at which scratch can be reused 49 Integer i_alpha ! scr ptr for exponents and Hermite expansion coef prefactor 50 Integer i_E ! Hermite expansion coefficients 51 Integer i_G1 ! center G coordinates (first use;through making of E) 52 Integer i_G2 ! center G coordinates (second use; making of R) 53 Integer i_Gt ! Temp array for G coorinates 54 Integer i_ABC2I ! Temp array for prefactors of E 55 Integer i_R0 ! the auxillary hermite integrals required for sum 56 Integer i_ijk ! pointer for auxillary hermite integrals 57 Integer i_RS ! prefactors for auxillary hermite integrals 58 Integer i_GC ! center for auxillary hermite integral recursion 59 Integer i_r ! the auxillary hermite integrals required for recursion 60 Integer i_ff ! prefactor array for forming recusion set of R's 61 62#if defined(INTDEBUG) 63 write(luout,*)' hf3pot debug' 64 if (.not.dryrun) then 65 write(luout,*)'b3nai_sz',b3nai_sz 66 write(luout,*)'dryrun',dryrun 67 write(luout,*)'lscr',lscr 68 call hf_print_set(1) 69 call hf_print('hf3pot a',Axyz,Aprims,Acoef,NPA,NCA,La) 70 call hf_print('hf3pot b',Bxyz,Bprims,Bcoef,NPB,NCB,Lb) 71 call hf_print('hf3pot c',Cxyz,Cprims,Ccoef,NPC,NCC,Lc) 72 call hf_print_set(0) 73 endif 74#endif 75c 76c Compute 3-center nuclear attraction integrals for ECPs 77c 78 gencon = nca.gt.1.or.ncb.gt.1.or.ncc.gt.1 79 if (gencon) call errquit 80 & (' hf3pot not ready for general contractions ',911, 81 & INT_ERR) 82c 83 i_alpha = 1 84 i_top = i_alpha + (NPA*NPB*NPC)*4 - 1 85 86#if defined(INTDEBUG) 87 write(luout,*)'i_alpha :',i_alpha 88 write(luout,*)'i_top :',i_top 89 write(luout,*)npa,npb,npc,(NPA*NPB*NPC),(NPA*NPB*NPC)*4 90#endif 91 if (i_top.gt.lscr) then 92 write(luout,*)' hf3pot: insufficient scratch space ' 93 write(luout,*)' : needed :',i_top 94 write(luout,*)' : allocated :',lscr 95 call errquit('hf3pot: fatal error ',1, INT_ERR) 96 endif 97c 98c.. determine actual product pairs to be kept. 99c The following line is to take care of compiler warnings. 100c 101 maxmem = i_top 102 if (DryRun) then 103 maxmem = i_top 104 Nabc = NPA*NPB*NPC 105 else 106 call hf1set3( 107 & Axyz,Aprims,Acoef,NPA, 108 & Bxyz,Bprims,Bcoef,NPB, 109 & Cxyz,Cprims,Ccoef,NPC, 110 & scr(i_alpha),Nabc) 111 endif 112#if defined(INTDEBUG) 113 write(luout,*) 114 write(luout,*)' nabc = ',nabc 115 write(luout,*) 116#endif 117 118c Define the center of the charge distribution. 119c A+B -> P; P + C -> G 120 121 Lg = La+Lb+Lc ! total angular momentum 122 123c 124 i_E = i_top + 1 125 i_G1 = i_E + Nabc*3*(Lg+1)*(La+1)*(Lb+1)*(Lc+1) 126 i_save = i_G1 ! i_alpha and I_E saved 127 i_top = i_G1 + Nabc*3 - 1 128 129#if defined(INTDEBUG) 130 write(luout,*) 131 write(luout,*)'i_alpha :',i_alpha 132 write(luout,*)'i_E :',i_E 133 write(luout,*)'i_G1 :',i_G1 134 write(luout,*)'i_save :',i_save 135 write(luout,*)'i_top :',i_top 136 write(luout,*) 137#endif 138 139 if (i_top.gt.lscr) then 140 write(luout,*)' hf3pot: insufficient scratch space ' 141 write(luout,*)' : needed :',i_top 142 write(luout,*)' : allocated :',lscr 143 144 write(luout,*)'i_alpha :',i_alpha 145 write(luout,*)'i_E :',i_E 146 write(luout,*)'i_G1 :',i_G1 147 call errquit('hf3pot: fatal error ',2, INT_ERR) 148 endif 149c 150 if (DryRun) then 151 maxmem = max(maxmem, i_top) 152 else 153 call hfctr3(Axyz,Bxyz,Cxyz,scr(i_alpha),scr(i_G1),Nabc) 154 endif 155c 156c Define the Hermite linear expansion coefficients 157c 158 i_Gt = i_top + 1 159 i_ABC2I = i_Gt + Nabc*3 160 i_top = i_ABC2I + Nabc*3 - 1 161 162#if defined(INTDEBUG) 163 write(luout,*) 164 write(luout,*)'i_save :',i_save 165 write(luout,*)'i_alpha :',i_alpha 166 write(luout,*)'i_E :',i_E 167 write(luout,*)'i_G1 :',i_G1 168 write(luout,*)'i_Gt :',i_Gt 169 write(luout,*)'i_ABC2I :',i_ABC2I 170 write(luout,*)'i_top :',i_top 171 write(luout,*) 172#endif 173 174 if (i_top.gt.lscr) then 175 write(luout,*)' hf3pot: insufficient scratch space ' 176 write(luout,*)' : needed :',i_top 177 write(luout,*)' : allocated :',lscr 178 179 write(luout,*)'i_alpha :',i_alpha 180 write(luout,*)'i_E :',i_E 181 write(luout,*)'i_G1 :',i_G1 182 write(luout,*)'i_Gt :',i_Gt 183 write(luout,*)'i_ABC2I :',i_ABC2I 184 call errquit('hf3pot: fatal error ',3, INT_ERR) 185 endif 186c 187 if (DryRun) then 188 maxmem = max(maxmem, i_top) 189 else 190 call hf1mke3(Axyz, Bxyz, Cxyz, 191 & scr(i_alpha), scr(i_G1), scr(i_Gt), scr(i_ABC2I), 192 & scr(i_E), Nabc, La, Lb, Lc) 193 endif 194c 195 Lg3 = ((Lg+1)*(Lg+2)*(Lg+3))/6 196 197#if defined(INTDEBUG) 198 write(luout,*) 199 write(luout,*)'lg,lg3',lg,lg3 200 write(luout,*) 201#endif 202 203 i_R0 = i_save 204 i_IJK = i_R0 + Nabc*Lg3 205 i_RS = i_IJK + (Lg+1)*(Lg+1)*(Lg+1) 206 i_GC = i_RS + Nabc 207 i_G2 = i_GC + Nabc*3 208 i_ff = i_G2 + Nabc*3 209 i_R = i_ff + 2*Nabc 210 i_top = i_R + Nabc*(Lg+1)*Lg3 - 1 211 212#if defined(INTDEBUG) 213 write(luout,*) 214 write(luout,*)'i_save :',i_save 215 write(luout,*)'i_alpha :',i_alpha 216 write(luout,*)'i_E :',i_E 217 write(luout,*)'i_R0 :',i_R0 218 write(luout,*)'i_IJK :',i_IJK 219 write(luout,*)'i_RS :',i_RS 220 write(luout,*)'i_GC :',i_GC 221 write(luout,*)'i_G2 :',i_G2 222 write(luout,*)'i_ff :',i_ff 223 write(luout,*)'i_R :',i_R 224 write(luout,*)'i_top :',i_top 225 write(luout,*) 226#endif 227 228 if (i_top.gt.lscr) then 229 write(luout,*)' hf3pot: insufficient scratch space ' 230 write(luout,*)' : needed :',i_top 231 write(luout,*)' : allocated :',lscr 232 233 write(luout,*)'i_alpha :',i_alpha 234 write(luout,*)'i_E :',i_E 235 write(luout,*)'i_R0 :',i_R0 236 write(luout,*)'i_IJK :',i_IJK 237 write(luout,*)'i_RS :',i_RS 238 write(luout,*)'i_GC :',i_GC 239 write(luout,*)'i_G2 :',i_G2 240 write(luout,*)'i_ff :',i_ff 241 write(luout,*)'i_R :',i_R 242 call errquit('hf3pot: fatal error ',3, INT_ERR) 243 endif 244c 245 if (DryRun) then 246 maxmem = max(maxmem, i_top) 247 else 248 call hf3mkr(Axyz, Bxyz, Cxyz, scr(i_alpha), scr(i_G2), 249 & scr(i_RS), scr(i_GC), scr(i_ff), scr(i_R), scr(i_R0), 250 & scr(i_IJK), Nabc, Lg, Lg3) 251 endif 252c 253 if (DryRun) then 254 lscr = maxmem 255 return 256 endif 257c 258 call hf3PEabc(b3nai, scr(i_E), scr(i_R0), scr(i_IJK), 259 & Nabc, La, Lb, Lc, Lg, Lg3, nint, b3nai_sz) 260c 261 end 262