1      Subroutine dim_grid_quadv0a(ncoef_max, natoms, iniz, ifin,
2     $                            rho, qwght, qxyz, xyz, expo, ccoef,
3     $                            Bmat, Fmat, Pmat, rq, cetobfr,
4     $                            ibf_ao, rqbuf, rchi_atom, g_vdim,
5     $                            nmat, do_scf, do_cpks_l, l3d,
6     $                            ipm, imag, dimxyz, muind)
7
8      implicit none
9#include "errquit.fh"
10#include "cgridfile.fh"
11#include "stdio.fh"
12#include "cdft.fh"
13#include "mafdecls.fh"
14#include "bas.fh"
15#include "global.fh"
16#include "util.fh"
17#include "steric.fh"
18#include "msgids.fh"
19#include "dimqm.fh"
20c
21c   Input Variables
22      double precision rqbuf(*)
23      double precision rad
24      integer g_vdim
25      double precision qxyz(3,n_per_rec)
26      double precision qwght(n_per_rec)
27      integer ncoef_max
28      integer natoms
29      integer iniz(natoms)
30      integer ifin(natoms)
31      double precision rho(*)
32      double precision xyz(3,natoms)
33      double precision expo(nbf_ao_mxprim)
34      double precision Bmat(nbf_ao_mxnbf_ce*n_per_rec)
35      double precision Pmat(*)
36      double precision Fmat(*)
37      double precision rq(n_per_rec, natoms)
38      integer cetobfr(2,natoms)
39      integer ibf_ao(nbf_ao)
40      double precision rchi_atom(natoms)
41      integer nmat
42      logical do_scf
43      logical do_cpks_l
44      double precision ccoef(ncoef_max)
45
46      integer iqsh
47      integer nqpts, ictr_buf
48      integer ncube,istep,ntot_cube,ncontrset
49      integer lbas_cent_info, ibas_cent_info,
50     &        lbas_cset_info, ibas_cset_info,
51     &        ldocset, idocset,i_iscratch,l_iscratch
52      integer ncontrsetx,lbas_cent_xinfo, ibas_cent_xinfo,
53     &        lbas_cset_xinfo, ibas_cset_xinfo,
54     .     ldocsetx, idocsetx
55
56      logical grid_file_rewind
57      external grid_file_rewind
58      logical xc_chkgrad, xc_chktau, kske
59      external xc_chkgrad, xc_chktau
60      logical l3d
61      integer ipm
62      integer imag
63      double precision dimxyz(3, nDIM)
64      double precision muind(3, nDIM, nmat)
65c
66c      if(ldebug) then
67c        write(luout,*) "dim_grid_quadv0a start"
68c      end if
69
70      if(n_rec_in_file.eq.0) goto 100
71c
72c     rewind grid pts file
73c
74      if (.not. grid_file_rewind())
75     $   call errquit('grid_quadv0a: rewinding gridpts?', 0,
76     &       UNKNOWN_ERR)
77      if (.not.bas_numcont(AO_bas_han, ncontrset))
78     &     call errquit('Exiting in grid_quadv0a',0, BASIS_ERR)
79c
80c     Allocate and create info for new basis function evaluator
81c
82      if (.not.MA_Push_Get(mt_int, 3*natoms, 'bas_cent_info',
83     &     lbas_cent_info, ibas_cent_info))
84     &     call errquit('grid_quadv0a: cannot allocate bas_cent_info',0,
85     &       MA_ERR)
86      if (.not.MA_Push_Get(mt_int, 6*ncontrset, 'bas_cset_info',
87     &     lbas_cset_info, ibas_cset_info))
88     &     call errquit('grid_quadv0a: cannot allocate bas_cset_info',0,
89     &       MA_ERR)
90
91      call xc_make_basis_info(AO_bas_han, int_mb(ibas_cent_info),
92     &     int_mb(ibas_cset_info), natoms)
93
94      if (.not.MA_Push_Get(mt_log, ncontrset, 'docset',
95     &     ldocset, idocset))
96     &     call errquit('grid_quadv0a: cannot allocate docset',
97     .     ncontrset, MA_ERR)
98c
99      ntot_cube=0
100c      write(luout,*) "n_recs_in_file:", n_rec_in_file
101      do 200 iqsh = 1, n_rec_in_file
102c      write(luout,*) "Top of loop", iqsh
103c
104c     Define the current range of radial shells and integration center.
105c
106         call grid_file_read(n_per_rec, nqpts, ictr_buf,
107     &        rad,rqbuf,nsubb)
108         if(nqpts.gt.buffer_size) call
109     '        errquit(' buffersize exceed by qpts ',nqpts, DISK_ERR)
110        if(nqpts.eq.0) goto 200
111        istep=0
112c        write(luout,*) iqsh, "nsubb:", nsubb
113        do  ncube=1,nsubb
114c
115c        temp put buf into currently used arrays qxyz and qwght
116c
117         call grid_repack(rqbuf, qxyz, qwght, nqpts, rad,istep)
118         if(nqpts.ne.0) then
119            call dim_grid_quadv0b(
120     $              ictr_buf, qxyz, qwght, nqpts, rad,
121     $              ncoef_max, natoms, iniz, ifin, rho,
122     $              xyz, expo, Bmat, Fmat, Pmat, rq,
123     $              cetobfr, ibf_ao, int_mb(ibas_cent_info),
124     $              int_mb(ibas_cset_info), log_mb(idocset), rchi_atom,
125     $              g_vdim, nmat, do_scf, do_cpks_l, l3d, ipm, imag,
126     $              dimxyz, muind)
127c
128             ntot_cube=ntot_cube+1
129c
130c
131         endif
132      enddo
133c      write(luout,*) "Bottom of loop", iqsh
134 200  continue
135      call ga_igop(Msg_Excrho, ntot_cube , 1, '+')
136c
137      if (.not.ma_chop_stack(lbas_cent_info))
138     &     call errquit('grid_quadv0a: pop stack failed.',1, MA_ERR)
139 100  continue
140c
141c      if(ldebug) then
142c        write(luout,*) "dim_grid_quadv0a end"
143c      end if
144      call ga_sync
145      return
146      end
147