1c
2c     Calculate the Hartree potential on the grid
3      subroutine gridHartreePotential(use_modelpotential,
4     &       ao_bas_han, geom, natoms, ipol, g_dens, npts, qxyz, qwght,
5     &       closegridpts, gexpo, gcoef, amat_coul)
6c
7      implicit none
8c
9#include "stdio.fh"
10c
11      integer npts,ao_bas_han,natoms,geom
12      integer g_dens(2),ipol
13      double precision qwght(npts)
14      double precision qxyz(3,npts)
15      integer i,j,k
16      double precision amat_coul(npts,ipol)
17      integer closegridpts(npts)
18      double precision tol
19      logical use_modelpotential
20      double precision gexpo(natoms,50)
21      double precision gcoef(natoms,50)
22c
23c     == use the model potential approach ==
24      if (use_modelpotential) then
25c
26c      == calculate the model density to construct the Coulomb contribution ==
27       call calc_modelpotential(geom,natoms,npts,qxyz,qwght,
28     &       closegridpts,gexpo,gcoef,amat_coul)  !amat_coul = v(modelpot)
29      else
30c
31c      == calculate the hartree potential on a supplied list of points ==
32       tol = 1d-8
33       do i=1,ipol
34          call potential_list(ao_bas_han, g_dens(i), npts, qxyz,
35     &                     amat_coul(1,i), tol)
36       enddo
37       do k = 1,npts
38          if (ipol.gt.1) amat_coul(k,1)=amat_coul(k,1)+amat_coul(k,2)
39          amat_coul(k,1) = -amat_coul(k,1)
40       enddo
41      end if ! use_modelpotential
42
43      return
44      end
45