1C 2C Initialization for charge density Coulomb fitting. 3C 4 subroutine rt_tddft_init_coulcdfit (params) 5 implicit none 6 7#include "errquit.fh" 8#include "mafdecls.fh" 9#include "stdio.fh" 10#include "global.fh" 11#include "msgids.fh" 12#include "util.fh" 13#include "cdft.fh" 14#include "rt_tddft.fh" 15#ifdef SCALAPACK 16 integer ga_cholesky,ga_llt_i 17 external ga_cholesky,ga_llt_i 18#endif 19 20 21C == In/out == 22 type(rt_params_t) params !cdfit params stored in here 23 24 25C == Parameters == 26 integer, parameter :: npol = 1 27 character(*), parameter :: pname = "rt_tddft_init_coulcdfit: " 28 29 30C == External == 31 logical dft_mem3c 32 external dft_mem3c 33 34 35C == Variables == 36 logical new_incore 37 integer n_batch 38 integer n3c_int, n3c_dbl, n_semi_bufs 39 integer iwhat_max 40 integer l_3ceri 41 integer k_3ceri 42 integer l_3cwhat 43 integer k_3cwhat 44 integer fd 45 46C (for part from dft_main0d.F) 47 integer g_tmpb 48 integer info 49 logical IOLGC 50 integer lmiss 51 integer me 52 53 54 me = ga_nodeid () 55 56 57C 58C (ripped from dft_main0d.F) 59C 60 IOLGC = .TRUE. 61 62C if (noio.eq.1) IOLGC = .FALSE. 63 64 if (params%nodisk) then 65 call errquit (pname//"not working with nodisk yet",0,0) 66 IOLGC = .false. 67 endif 68 69 70C call errquit (pname//"XXX check npol, always 1?",0,0) 71C call rt_tddft_print_warning ("CHECK THAT CD FITTING WORKING") 72 73 74 if (CDFIT)then 75c 76c Determine the characteristics of the CD Gaussian basis set. 77c 78c 79c Compute the matrix inverse of the CD 2-ctr ERIs. 80c 81 if (.not. ga_create(mt_dbl, nbf_cd, nbf_cd, 'CD 2cERI', 82 & 0, nbf_cd, g_2ceri)) 83 & call errquit(pname//'Error creating g_2ceri',0, 84 & GA_ERR) 85 call ga_zero(g_2ceri) 86 call dft_get2eri(CD_bas_han, g_2ceri,oskel) 87 if (oskel)call 88 . sym_symmetrize(geom,cd_bas_han,.false.,g_2ceri) 89 call ga_sync() 90 if (.not. ga_duplicate(g_2ceri, g_cdinv, 'CD 2cERInv')) 91 & call errquit(pname//'Error creating g_cdinv',0, GA_ERR) 92 call ga_zero(g_cdinv) 93 lmiss = 1 94c if (lmiss.eq.1) then 95 call ga_zero(g_cdinv) 96 info = 0 97#if defined(PARALLEL_DIAG) 98#ifdef SCALAPACK 99 call ga_copy(g_2ceri, g_cdinv) 100 call ga_sync() 101 info= ga_cholesky('U',g_cdinv) 102#else 103 call ga_chol(g_2ceri, g_cdinv, info) 104#endif 105#else 106 call ga_chol_seq(g_2ceri, g_cdinv, info) 107#endif 108 if (info.ne.0)then 109 if (me.eq.0)then 110 write(LuOut,*)' Problem in performing a Cholesky ' 111 write(LuOut,*)' decomposition of the 2-ctr ERI ' 112 write(LuOut,*)' matrix using CD fitting basis. ' 113 write(LuOut,*)' Attempting a diag/inverse. ' 114 endif 115 endif 116 if (info.eq.0) then 117 g_tmpb = g_2ceri 118#if defined(PARALLEL_DIAG) 119#ifdef SCALAPACK 120 info = ga_llt_i('U',g_cdinv,-1) 121 if (info.ne.0)then 122 if (me.eq.0)then 123 write(LuOut,*)' Problem in performing a Invers. ' 124 write(LuOut,*)' of the 2-ctr ERI ' 125 endif 126 call ga_sync 127 call errquit(pname//'Inverse failed ',0,0) 128 endif 129 130#else 131 call ga_inverse(g_cdinv, g_tmpb) 132#endif 133#ifndef SCALAPACK 134 call ga_dgemm('T', 'N', nbf_cd, nbf_cd, nbf_cd, 1.d0, 135 & g_tmpb, g_tmpb, 0.d0, g_cdinv) 136#endif 137#else 138 call ga_copy(g_cdinv, g_tmpb) 139 call ga_inv_seq(g_tmpb, g_cdinv) 140#endif 141 else 142 call dft_invdiag(g_2ceri, g_cdinv, 143 & nbf_cd) 144 endif 145#ifndef SCALAPACK 146c 147c second build of g_2ceri needed becuase previous calls destroyed it 148c 149 call ga_zero(g_2ceri) 150 call dft_get2eri(CD_bas_han, g_2ceri,oskel) 151 if (oskel)call 152 . sym_symmetrize(geom,cd_bas_han,.false.,g_2ceri) 153#endif 154 if (IOLGC.and.(me.eq.0)) then 155 lmiss = 0 156 call dft_invio('CDI', g_cdinv, nbf_cd, 'WRITE', lmiss) 157 if (lmiss.ne.0)call errquit 158 & (pname//'dft_invio - abnormal write of CDI ', 0, 159 & DISK_ERR) 160 lmiss = 0 161 call dft_invio('CD', g_2ceri, nbf_cd, 'WRITE', lmiss) 162 if (lmiss.ne.0)call errquit 163 & (pname//'dft_invio - abnormal write of CD ', 0, 164 & DISK_ERR) 165 endif 166c endif 167 if (IOLGC) then 168 if (.not. ga_destroy(g_cdinv)) call errquit 169 & (pname//'Could not destroy g_xcinv', 0, GA_ERR) 170 if (.not. ga_destroy(g_2ceri)) call errquit 171 & (pname//'Could not destroy g_xcinv', 0, GA_ERR) 172 endif 173 endif 174C 175C (end from dft_main0d.F) 176C 177 178 179C 180C Set up three center integrals. This was already done in SCF, but 181C we have to redo it here. 182C 183 if (dft_mem3c(params%rtdb, 184 $ params%natoms, npol, .false., .false., 185 $ n3c_int, n3c_dbl, !! n_semi_bufs, 186 $ l_3ceri, k_3ceri, l_3cwhat, k_3cwhat)) then 187 188 incore=.false. 189 call dft_3cincor (n_batch, n3c_int, int_mb(k_3cwhat), 190 $ dbl_mb(k_3cERI), n3c_dbl) 191 192 new_incore = .true. 193 else 194 new_incore = .false. 195 endif 196 197 198C 199C Replace the "incore" flag in cdft.fh 200C 201 incore = new_incore 202 203C if (new_incore .neqv. incore) 204C $ call errquit (pname//"incore different from SCF", 0, 0) 205 206 207C 208C Load params necessary for CD fitting into "params" struct. 209C 210 params%n_batch = n_batch 211 params%k_3ceri = k_3ceri 212 params%l_3ceri = l_3ceri 213 params%n3c_int = n3c_int 214 params%k_3cwhat = k_3cwhat 215 params%l_3cwhat = l_3cwhat 216 params%n3c_dbl = n3c_dbl 217 params%iwhat_max = iwhat_max 218 params%n_semi_bufs = n_semi_bufs 219 params%fd = fd 220 221 222 end subroutine rt_tddft_init_coulcdfit 223 224 225 226c $Id$ 227