1!! gen1int: compute one-electron integrals using rotational London atomic-orbitals 2!! Copyright 2009-2012 Bin Gao, and Andreas Thorvaldsen 3!! 4!! gen1int is free software: you can redistribute it and/or modify 5!! it under the terms of the GNU Lesser General Public License as published by 6!! the Free Software Foundation, either version 3 of the License, or 7!! (at your option) any later version. 8!! 9!! gen1int is distributed in the hope that it will be useful, 10!! but WITHOUT ANY WARRANTY; without even the implied warranty of 11!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 12!! GNU Lesser General Public License for more details. 13!! 14!! You should have received a copy of the GNU Lesser General Public License 15!! along with gen1int. If not, see <http://www.gnu.org/licenses/>. 16!! 17!! Tests subroutines in contr_cgto_carmom.F90 18!! 19!! 2011-04-03, Bin Gao 20!! * first version 21 22 !> \brief tests subroutines in contr_cgto_carmom.F90 23 !> \author Bin Gao 24 !> \date 2011-04-03 25 !> \param io_log is the IO unit of log file 26 !> \return test_failed indicates if the test is failed 27 subroutine test_contr_cgto_carmom(io_log, test_failed) 28 use xkind 29 ! module of HTML test log routines 30 use html_log 31 implicit none 32 integer, intent(in) :: io_log 33 logical, intent(out) :: test_failed 34! water molecule using cc-pV6Z basis sets 35#include "test_f90/water_ccpv6z.h" 36! test parameters and referenced results 37#include "water_ccpv6z_cgto_carmom.h" 38 integer itst, jtst, ktst, ltst !incremental recorders over test cases 39 integer iblock, jblock !incremental recorders over blocks of sub-shells 40 integer num_opt !number of operators including derivatives 41 real(REALK), allocatable :: contr_ints(:,:,:,:,:) 42 !contracted integrals 43 integer idx_bra, idx_ket !indices of bra and ket centers 44 integer num_cgto_bra, num_cgto_ket !number of Cartesian Gaussians of bra and ket centers 45 integer ierr !error information 46 integer addr_ref !address of referenced results 47 real(REALK) begin_time !begin of CPU time 48 real(REALK) end_time !end of CPU time 49 logical different !if result from \fn(contr_cgto_carmom) is different from reference 50 real(REALK), allocatable :: ro_cints(:,:,:,:,:) !reordered integrals 51 integer, allocatable :: power_bra(:,:) !powers on bra and ket centers 52 integer, allocatable :: power_ket(:,:) 53 integer xpow, ypow !incremental recorders over xy powers 54 integer addr_bra, addr_ket !addresses of GTOs on bra and ket centers 55 integer icontr, jcontr, icgto, jcgto, iopt !incremental recorders 56 ! we first test subroutines in reorder_ints.F90 which will be used later on 57 call html_log_text("We first test subroutines in reorder_ints.F90 ...", io_log) 58 call test_reorder_ints(io_log, test_failed) 59 if (test_failed) then 60 call html_log_text(log_text="Tests subroutines in reorder_ints.F90 failed! "// & 61 "No further test will be performed!", & 62 io_log=io_log, font_color="red") 63 return 64 else 65 call html_log_text("Tests subroutines in reorder_ints.F90 passed!", io_log) 66 end if 67 ! dumps information of tests 68 call html_log_text("Tests water molecule using cc-pV6Z Cartesian GTOs ...", io_log) 69 ! loops over test cases 70 jtst = 1 71 addr_ref = 0 72 do itst = 1, NUM_TESTS 73 ktst = jtst+2 74 call html_log_int_number("Order of electronic derivatives:", & 75 ORDER_ELEC(itst), "I4", io_log) 76 if (IDX_DIPORG(itst)>0) & 77 call html_log_int_number("Index of dipole origin:", & 78 IDX_DIPORG(itst), "I4", io_log) 79 call html_log_real_array("Coordinates of dipole origin:", DIP_ORIGIN(jtst:ktst), & 80 "F16.8", io_log) 81 call html_log_real_number("Scaling constant:", SCAL_CONST(itst), "F16.8", io_log) 82 call html_log_int_number("Order of Cartesian multipole moments:", & 83 ORDER_MOM(itst), "I4", io_log) 84 call html_log_int_number("Order of partial geometric derivatives on bra center:", & 85 ORDER_GEO_BRA(itst), "I4", io_log) 86 call html_log_int_number("Order of partial geometric derivatives on ket center:", & 87 ORDER_GEO_KET(itst), "I4", io_log) 88 call html_log_int_number("Order of geometric derivatives on dipole origin:", & 89 ORDER_GEO_MOM(itst), "I4", io_log) 90 ! computes the number of total geometric derivatives 91 num_opt = 1 92 do ltst = jtst, jtst+NUM_CENTS(itst)-1 93 num_opt = num_opt*(ORDER_CENT(ltst)+1)*(ORDER_CENT(ltst)+2)/2 94 end do 95 ltst = jtst+NUM_CENTS(itst)-1 96 if (NUM_CENTS(itst)>0) then 97 call html_log_int_number("Numbers of differentiated centers of total "// & 98 "geometric derivatives:", NUM_CENTS(itst), "I4", io_log) 99 call html_log_int_array("Indices of differentiated centers:", & 100 IDX_CENT(jtst:ltst), "I4", io_log) 101 call html_log_int_array("Orders of differentiated centers:", & 102 ORDER_CENT(jtst:ltst), "I4", io_log) 103 if (NUM_CENTS(itst)>3) then 104 call html_log_int_number(log_text="Invalid number of differentiated "// & 105 "centers of total geometric derivatives:", & 106 log_int=NUM_CENTS(itst), fmt_int="I4", & 107 io_log=io_log, font_color="red") 108 test_failed = .true. 109 cycle 110 end if 111 end if 112 ! gets the begin time 113 call xtimer_set(begin_time) 114 ! computs the number of operators including the derivatives 115 num_opt = num_opt & 116 * (ORDER_ELEC(itst)+1)*(ORDER_ELEC(itst)+2) & 117 * (ORDER_MOM(itst)+1)*(ORDER_MOM(itst)+2) & 118 * (ORDER_GEO_BRA(itst)+1)*(ORDER_GEO_BRA(itst)+2) & 119 * (ORDER_GEO_KET(itst)+1)*(ORDER_GEO_KET(itst)+2) & 120 * (ORDER_GEO_MOM(itst)+1)*(ORDER_GEO_MOM(itst)+2)/32 121 ! gets the blocks of sub-shells 122 iblock = BRA_BLOCK(itst) 123 jblock = KET_BLOCK(itst) 124 call html_log_int_number("ID of block of sub-shells on bra center:", & 125 iblock, "I4", io_log) 126 call html_log_int_number("ID of block of sub-shells on ket center:", & 127 jblock, "I4", io_log) 128 ! gets the indices of bra and ket centers 129 idx_bra = IDX_CENT_BLOCK(iblock) 130 idx_ket = IDX_CENT_BLOCK(jblock) 131 ! sets the number of Cartesian Gaussians of bra and ket centers 132 num_cgto_bra = (ANGULAR_BLOCK(iblock)+1)*(ANGULAR_BLOCK(iblock)+2)/2 133 num_cgto_ket = (ANGULAR_BLOCK(jblock)+1)*(ANGULAR_BLOCK(jblock)+2)/2 134 ! allocates the contracted integrals 135 allocate(contr_ints(num_cgto_bra,NUM_CONTR_BLOCK(iblock), & 136 num_cgto_ket,NUM_CONTR_BLOCK(jblock), & 137 num_opt), stat=ierr) 138 if (ierr/=0) then 139 call html_log_int_number(log_text="Failed to allocate contr_ints:", & 140 log_int=num_cgto_bra*NUM_CONTR_BLOCK(iblock) & 141 *num_cgto_ket*NUM_CONTR_BLOCK(jblock) & 142 *num_opt, & 143 fmt_int="I12", io_log=io_log, font_color="red") 144 test_failed = .true. 145 cycle 146 end if 147 ! computes the contracted integrals 148 call contr_cgto_carmom(idx_bra, COORDS(3*idx_bra-2:3*idx_bra), & 149 ANGULAR_BLOCK(iblock), NUM_PRIM_BLOCK(iblock), & 150 EXPONENT_BLOCK(START_PRIM_BLOCK(iblock)+1 & 151 :START_PRIM_BLOCK(iblock+1)), & 152 NUM_CONTR_BLOCK(iblock), & 153 REF_NRM_CGTO(START_CONTR_BLOCK(iblock)+1 & 154 :START_CONTR_BLOCK(iblock+1)), & 155 idx_ket, COORDS(3*idx_ket-2:3*idx_ket), & 156 ANGULAR_BLOCK(jblock), NUM_PRIM_BLOCK(jblock), & 157 EXPONENT_BLOCK(START_PRIM_BLOCK(jblock)+1 & 158 :START_PRIM_BLOCK(jblock+1)), & 159 NUM_CONTR_BLOCK(jblock), & 160 REF_NRM_CGTO(START_CONTR_BLOCK(jblock)+1 & 161 :START_CONTR_BLOCK(jblock+1)), & 162 ORDER_ELEC(itst), IDX_DIPORG(itst), & 163 DIP_ORIGIN(jtst:ktst), SCAL_CONST(itst), & 164 ORDER_MOM(itst), ORDER_GEO_BRA(itst), & 165 ORDER_GEO_KET(itst), ORDER_GEO_MOM(itst), & 166 NUM_CENTS(itst), IDX_CENT(jtst:ltst), & 167 ORDER_CENT(jtst:ltst), num_cgto_bra, & 168 num_cgto_ket, num_opt, contr_ints) 169 ! gets the end time 170 call xtimer_set(end_time) 171 ! prints the CPU elapsed time 172 call html_log_real_number(log_text="Time (s) used for contr_cgto_carmom:", & 173 log_real=end_time-begin_time, fmt_real="F10.4", & 174 io_log=io_log, font_color="blue") 175 ! prepares the powers on bra and ket centers 176 allocate(power_bra(3,num_cgto_bra), stat=ierr) 177 if (ierr/=0) then 178 call html_log_int_number(log_text="Failed to allocate power_bra:", & 179 log_int=3*num_cgto_bra, fmt_int="I12", & 180 io_log=io_log, font_color="red") 181 deallocate(contr_ints) 182 test_failed = .true. 183 cycle 184 end if 185 allocate(power_ket(3,num_cgto_ket), stat=ierr) 186 if (ierr/=0) then 187 call html_log_int_number(log_text="Failed to allocate power_ket:", & 188 log_int=3*num_cgto_ket, fmt_int="I12", & 189 io_log=io_log, font_color="red") 190 deallocate(power_bra) 191 deallocate(contr_ints) 192 test_failed = .true. 193 cycle 194 end if 195 addr_bra = 0 196 do xpow = ANGULAR_BLOCK(iblock), 0, -1 197 do ypow = ANGULAR_BLOCK(iblock)-xpow, 0, -1 198 addr_bra = addr_bra+1 199 power_bra(1,addr_bra) = xpow 200 power_bra(2,addr_bra) = ypow 201 power_bra(3,addr_bra) = ANGULAR_BLOCK(iblock)-(xpow+ypow) 202 end do 203 end do 204 addr_ket = 0 205 do xpow = ANGULAR_BLOCK(jblock), 0, -1 206 do ypow = ANGULAR_BLOCK(jblock)-xpow, 0, -1 207 addr_ket = addr_ket+1 208 power_ket(1,addr_ket) = xpow 209 power_ket(2,addr_ket) = ypow 210 power_ket(3,addr_ket) = ANGULAR_BLOCK(jblock)-(xpow+ypow) 211 end do 212 end do 213 ! reorders the integrals in the order of Dalton format 214 allocate(ro_cints(num_cgto_bra,NUM_CONTR_BLOCK(iblock), & 215 num_cgto_ket,NUM_CONTR_BLOCK(jblock), & 216 num_opt), stat=ierr) 217 if (ierr/=0) then 218 call html_log_int_number(log_text="Failed to allocate ro_cints:", & 219 log_int=num_cgto_bra*NUM_CONTR_BLOCK(iblock) & 220 *num_cgto_ket*NUM_CONTR_BLOCK(jblock) & 221 *num_opt, & 222 fmt_int="I12", io_log=io_log, font_color="red") 223 deallocate(power_bra) 224 deallocate(power_ket) 225 deallocate(contr_ints) 226 test_failed = .true. 227 cycle 228 end if 229 call reorder_cgto_ints(ANGULAR_BLOCK(iblock), num_cgto_bra, power_bra, & 230 ANGULAR_BLOCK(jblock), num_cgto_ket, power_ket, & 231 NUM_CONTR_BLOCK(iblock), NUM_CONTR_BLOCK(jblock), & 232 num_opt, contr_ints, ro_cints) 233 ! cleans 234 deallocate(contr_ints) 235 deallocate(power_bra) 236 deallocate(power_ket) 237 ! checks the results 238 do iopt = 1, num_opt 239 do jcontr = 1, NUM_CONTR_BLOCK(jblock) 240 do jcgto = 1, num_cgto_ket 241 do icontr = 1, NUM_CONTR_BLOCK(iblock) 242 do icgto = 1, num_cgto_bra 243 addr_ref = addr_ref+1 244 call check_difference(REF_CONTR_INTS(addr_ref), & 245 ro_cints(icgto,icontr,jcgto,jcontr,iopt), & 246 different) 247 if (different) then 248 call html_log_int_number("ID of test:", itst, "I3", io_log) 249 call html_log_int_array(log_text="ID of contraction and operator:", & 250 log_int=(/icgto,icontr,jcgto,jcontr,iopt/), & 251 fmt_int="I3", io_log=io_log) 252 call html_log_real_number( & 253 log_text="Referenced contracted integral:", & 254 log_real=REF_CONTR_INTS(addr_ref), & 255 fmt_real="Es16.8", io_log=io_log, font_color="blue") 256 call html_log_real_number( & 257 log_text="Result from contr_cgto_carmom:", & 258 log_real=ro_cints(icgto,icontr,jcgto,jcontr,iopt), & 259 fmt_real="Es16.8", io_log=io_log, font_color="red") 260 test_failed = .true. 261 end if 262 end do 263 end do 264 end do 265 end do 266 end do 267 ! cleans 268 deallocate(ro_cints) 269 ! updates the incremental recorder 270 jtst = jtst+3 271 end do 272 return 273 end subroutine test_contr_cgto_carmom 274