1 subroutine esp_esp(ene,occ,dns,x,q,grid,val) 2c 3c $Id$ 4c 5 implicit none 6c 7#include "esp_params.fh" 8#include "esp_common.fh" 9#include "global.fh" 10#include "bas.fh" 11#include "rtdb.fh" 12#include "geom.fh" 13#include "mafdecls.fh" 14#include "msgids.fh" 15#include "util.fh" 16#include "schwarz.fh" 17#include "stdio.fh" 18c 19 integer ga_create_atom_blocked 20 external ga_create_atom_blocked 21c 22 real*8 ene(2*nbf),occ(2*nbf),dns(mbf,mbf) 23 real*8 x(3,natoms),q(natoms) 24 real*8 grid(3,mxgrid),val(mxgrid) 25c 26 integer i,j,ish,jsh,ix 27 integer ilo,ihi,jlo,jhi 28 real*8 dist,vt 29 integer imin,nval,jshi 30 real*8 dmin,fact 31 character*10 today,now 32c 33 integer maxbuf,maxscr,avail 34 integer batch,no_batches,ngrid_batches 35 integer iptr 36 logical ldummy(1) 37c 38c get the density matrix, occupation and energies 39c 40 call esp_denmat(occ,ene) 41c 42 if (lscrn.ne.0) call schwarz_init(igeom,ibasis) 43c 44c get electrostatic potential on the grid points 45c 46 do i=1,ngrid 47 val(i)=0.0d0 48 end do 49c 50c estimate mem. req. 51c 52 call int_init_1eelec(maxbuf,maxscr,ibasis,0,ngrid) 53 avail=ma_inquire_avail(mt_dbl)*8/10 54 no_batches=nint(dble(maxbuf+maxscr)/dble(avail)) 55 no_batches=max(no_batches,1) 56 ngrid_batches=ngrid/no_batches 57#ifdef DEBUG 58 if(me.eq.0) then 59 write(6,*) ' avail ',avail 60 write(6,*) ' no_batches ',no_batches 61 write(6,*) ' ngrid_batches ',ngrid_batches 62 endif 63#endif 64 iptr=1 65 do batch=1,no_batches 66 if(batch.eq.no_batches) ngrid_batches=ngrid-iptr+1 67 call hnd_elfcon_schw(ibasis,igeom,lg_d, 68 C grid(1,iptr),ngrid_batches,val(iptr),0) 69 iptr=iptr+ngrid_batches 70 enddo 71 if(me.eq.0) then 72 do 8 i=1,ngrid 73 imin=0 74 val(i)=-val(i) 75 vt=val(i) 76 do 9 ix=1,natoms 77 dist=sqrt((grid(1,i)-x(1,ix))*(grid(1,i)-x(1,ix))+ 78 + (grid(2,i)-x(2,ix))*(grid(2,i)-x(2,ix))+ 79 + (grid(3,i)-x(3,ix))*(grid(3,i)-x(3,ix))) 80 if(imin.eq.0) then 81 dmin=dist 82 imin=ix 83 endif 84 if(dmin.gt.dist) then 85 dmin=dist 86 imin=ix 87 endif 88 val(i)=val(i)+q(ix)/dist 89 9 continue 90 8 continue 91 endif 92c 93 if(lscrn.ne.0) call schwarz_tidy() 94c 95 if(np.gt.1) then 96 call ga_brdcst(mre_002,val,ngrid*ma_sizeof(mt_dbl,1,mt_byte),0) 97 endif 98c 99c deallocate memory for the integrals 100c 101 if(.not.ga_destroy(lg_d)) 102 + call md_abort('esp_denmat: ga_destroy lg_d failed',me) 103c 104c integral termination 105c 106 call int_terminate() 107c 108 if(me.eq.0) then 109 call util_file_name('grid',.false.,.false.,grdfil) 110 open(unit=lfngrd,file=grdfil,form='formatted',status='unknown') 111 rewind(lfngrd) 112 write(lfngrd,'(i10,f20.10)') ngrid,charge 113 do i=1,ngrid 114 write(lfngrd,'(4f20.10)') (grid(j,i),j=1,3),val(i) 115 enddo 116 close(unit=lfngrd) 117 endif 118c 119 return 120 end 121