1 Subroutine nbf_to_mbf(basis_hand, mbf, mcenters, 2 & ctr_quad, rad_quad, xyz, 3 & ibf,iniz, ifin, zprim, acc_gauss) 4c 5C$Id$ 6c 7 implicit none 8#include "errquit.fh" 9c 10 integer basis_hand 11 integer mcenters ! [input] 12 double precision acc_gauss ! [input] 13 double precision xyz(3,*) ! [input] 14 integer ctr_quad ! grid center [input] 15 integer mbf ! [input/output] 16 integer ibf(*) ! [output] 17 double precision rad_quad ! [input] 18c 19#include "bas.fh" 20c 21c Distance Squared between Sampling Points and Centers 22c 23 integer iniz(mcenters),ifin(mcenters) 24 double precision zprim(*) 25 integer ncontrset, n1, icset, ictr, nprim, ncontr, 26 & isphere, l, iprimo,npt,nshbf 27 double precision zmin,acc_loc,r2,x,y,z 28 integer n,ic1,ic2,m 29 double precision alpha,logeps 30 double precision gaussian_range,r_arg 31 gaussian_range(n,logeps,alpha) = 32 $ (n*log(-logeps) - n*log(alpha) - 4.0d0*logeps) / 33 $ sqrt(-16.0d0*alpha*logeps) 34c 35 if (.not.bas_numcont(basis_hand, ncontrset)) 36 & call errquit('Exiting in nbf_to_mbf.',1, BASIS_ERR) 37c 38 call ifill(mcenters,0,iniz,1) 39 call ifill(mcenters,0,ifin,1) 40 acc_loc=-acc_gauss 41c 42 n1 = 0 43 npt=0 44 do 400 ictr=1,mcenters 45 if (.not.bas_ce2cnr(basis_hand,ictr,ic1,ic2)) 46 & call errquit('Exiting in xc_signf.',11, BASIS_ERR) 47 r2=0d0 48 if (ictr.ne.ctr_quad) then 49 x = xyz(1,ctr_quad) - xyz(1,ictr) 50 y = xyz(2,ctr_quad) - xyz(2,ictr) 51 z = xyz(3,ctr_quad) - xyz(3,ictr) 52 r2 = sqrt(x*x + y*y + z*z) 53 endif 54 r_arg=0d0 55 if (rad_quad.lt.r2) r_arg = (r2-rad_quad) 56 do icset = ic1,ic2 57c 58c get info about current contraction set 59c 60 if (.not. bas_continfo(basis_hand, icset, l ,nprim, 61 & ncontr, isphere)) 62 & call errquit('Exiting in xc_signf.',4, BASIS_ERR) 63c 64c 65c get exponents and contraction coefficients for this contraction set 66c 67 if (.not.bas_get_exponent(basis_hand, icset, zprim)) 68 & call errquit('Exiting in xc_signf.',7, BASIS_ERR) 69c 70c Determine the minimum Gaussian exponent. 71c 72 73 zmin = 1.D+06 74 do iprimo = 1,nprim 75 zmin = min(zprim(iprimo),zmin) 76 enddo 77c 78c Only include those basis functions that are "non-zero" for at least 79c one point in the sampling set. 80c 81 if(isphere.eq.1) then 82 nshbf=ncontr*(2*l+1) 83 else 84 nshbf=ncontr*(((l+1)*(l+2))/2) 85 endif 86 if (r_arg.lt.gaussian_range(l,acc_loc,zmin)) then 87 do m=1,nshbf 88 ibf(n1+m ) = npt+m 89 enddo 90 if (iniz(ictr).eq.0) iniz(ictr)=n1+1 91 n1=n1+nshbf 92 endif 93 94 npt=npt+nshbf 95 enddo 96c 97 ifin(ictr)= n1 98c 99 400 continue 100c 101 mbf = n1 102c 103 return 104 end 105