1! 2! Dalton, a molecular electronic structure program 3! Copyright (C) by the authors of Dalton. 4! 5! This program is free software; you can redistribute it and/or 6! modify it under the terms of the GNU Lesser General Public 7! License version 2.1 as published by the Free Software Foundation. 8! 9! This program 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 GNU 12! Lesser General Public License for more details. 13! 14! If a copy of the GNU LGPL v2.1 was not distributed with this 15! code, you can obtain one at https://www.gnu.org/licenses/old-licenses/lgpl-2.1.en.html. 16! 17! 18module mlcc3_active_spaces 19! 20! 21! mlcc3 initialiasion of active space variables 22! Authors Henrik Koch and Rolf H. Myhre 23! December 2014 24! 25 use mlcc_typedef 26 use mlcc_block_import 27 use mlcc_work 28 use mlcc3_data 29! 30! 31 implicit none 32! 33! Active space variables 34! 35 integer :: n_spaces, cc3_spaces, include_spaces 36 integer :: n_occ_act, n_vir_act, n_occ_gen, n_vir_gen 37 integer :: n_v_a_3 38 integer :: n_occ_inact, n_vir_inact, n_gen_inact 39 integer :: n_vir_aag, n_occ_aag !active*active*general 40 integer :: n_occ_integral, n_vir_integral !number of integrals 41 integer, dimension(:), pointer :: n_occ_space, n_vir_space 42! 43! CVS variables 44! 45 logical :: cvs_log, rmc_log, ion_log 46 integer :: n_cvs_rmc, n_ion 47 integer, dimension(:), pointer :: cvs_rmc_index, ion_index 48! 49contains 50! 51 subroutine center_h_import 52! 53 implicit none 54! 55 integer, parameter :: mxbsce = 2000, mxacat = 50, mxacbs = 5000 56 integer, parameter :: mxexbs = 500 57 integer, parameter :: mxspa = 5, maxspa = mxspa+1 58 integer, parameter :: mxcent = 500 59! 60 integer :: allocate_ok, i 61! 62 double precision thacoc, thacvi 63! 64! 65 integer nbcent, ibcent, nbscen, kbscen, & 66 & nacatm, lactat, lacbas, nactbs, nabsto, & 67 & nactoc, nactvi, ninaoc, ninavi, & 68 & ioract, noract, nactfr, nextbs, iextbs, & 69 & nocvec, nvivec, nacinp, lacinp, iordec, & 70 & nseld, iseld, nfract, & 71 & nabsoc, nabsvi, nacbsv,lacbsv, & 72 & nspace, natoac, labspa, nspaoc, nspavi, & 73 & nlfthf, nlftvi, iof2hf, iof2vi, iepshf, iepsvi, & 74 & nabso2, lacba2, nabsv2, lacbv2, iexpvi, iexpoc 75! 76 logical actsel, atomic, actfre, difadd, nboexp, seldir, & 77 & opnshl, fuldec, dospread, minspr, & 78 & dialst, extern, newact, spdils, loconl, addorb, & 79 & addexp 80! 81 common /center/ thacoc, thacvi, & 82 & nbcent(mxcent,8), ibcent(mxcent,8), & 83 & nbscen(mxcent), kbscen(mxbsce,mxcent), & 84 & nacatm, lactat(mxacat), lacbas(mxacbs), & 85 & nactbs(8), nabsto, nactfr, & 86 & ioract(8), noract(8), & 87 & nactoc(8), nactvi(8), ninaoc(8), ninavi(8), & 88 & nextbs(8), iextbs(mxexbs,8), & 89 & nocvec(mxcent,8), nvivec(mxcent,8), & 90 & nacinp, lacinp(mxacat), iordec(mxcent), & 91 & nseld(mxacat), iseld(mxbsce,mxacat), nfract(8), & 92 & nspace, natoac(mxspa), labspa(mxacat,mxspa), & 93 & nspaoc(8,maxspa), nspavi(8,maxspa), & 94 & nlfthf(8), nlftvi(8), iof2hf(8,maxspa), & 95 & iof2vi(8,maxspa), iepshf(8,maxspa), & 96 & iepsvi(8,maxspa), & 97 & nabsoc, nabsvi, nacbsv(8),lacbsv(mxacbs), & 98 & nabso2(mxspa), lacba2(mxacbs,mxspa), & 99 & nabsv2(mxspa), lacbv2(mxacbs,mxspa), & 100 & actsel, atomic, actfre, difadd, nboexp, seldir, & 101 & opnshl, fuldec, dospread, minspr, & 102 & dialst, extern, newact, spdils, loconl, & 103 & addorb, addexp, iexpvi, iexpoc 104! 105! 106! CVS, RMC, and ionisation not implemented in this version 107! Set to to false 108 cvs_log = .false. 109 ion_log = .false. 110 rmc_log = .false. 111! 112 if (mlcc3_active .and. .not. mlcc3_nrg_spa) then 113! 114 n_spaces = nspace + 1 115! 116 write(lupri,*) 117 write(lupri,*) 'nspace', nspace 118 write(lupri,*) 'n_spaces', n_spaces 119 write(lupri,*) 120! 121 if(n_general .gt. n_spaces) then 122! 123 call quit('n_general greater than n_spaces in mlcc3 center import') 124! 125 end if 126! 127! 128 allocate(n_occ_space(n_spaces), stat = allocate_ok) 129 if(allocate_ok /= 0) call quit('Something wrong with allocation in center_h_import') 130! 131 allocate(n_vir_space(n_spaces), stat = allocate_ok) 132 if(allocate_ok /= 0) call quit('Something wrong with allocation in center_h_import') 133! 134 do i = 1,n_spaces 135 n_occ_space(i) = nspaoc(1,i) 136 n_vir_space(i) = nspavi(1,i) 137 enddo 138! 139 n_occ_act = 0 140 n_vir_act = 0 141 n_occ_gen = 0 142 n_vir_gen = 0 143! 144 cc3_spaces = n_active ! Spaces treated with CC3 or higher 145 include_spaces = n_general ! Spaces to include for general orbital indexes 146! 147 do i = 1,cc3_spaces 148 n_occ_act = n_occ_act + n_occ_space(i) 149 n_vir_act = n_vir_act + n_vir_space(i) 150 end do 151 do i = 1,include_spaces 152 n_occ_gen = n_occ_gen + n_occ_space(i) 153 n_vir_gen = n_vir_gen + n_vir_space(i) 154 end do 155! 156 n_occ_inact = n_occ - n_occ_act 157 n_vir_inact = n_vir - n_vir_act 158! 159 n_gen_inact = n_occ - n_occ_gen 160! 161 else if(mlcc3_nrg_spa) then 162! 163 n_occ_act = n_occ_inp 164 n_vir_act = n_vir_inp 165! 166 if(mlcc3_nrg_gen) then 167 n_occ_gen = n_gen_o_inp 168 n_vir_gen = n_gen_v_inp 169 else 170 n_occ_gen = n_occ 171 n_vir_gen = n_vir 172 end if 173! 174 n_occ_inact = n_occ - n_occ_act 175 n_vir_inact = n_vir - n_vir_act 176! 177 n_gen_inact = n_occ - n_occ_gen 178! 179 else 180! 181 n_occ_act = n_occ 182 n_occ_gen = n_occ 183 n_vir_act = n_vir 184 n_vir_gen = n_vir 185! 186 n_occ_inact = 0 187 n_vir_inact = 0 188! 189 n_gen_inact = 0 190! 191 end if 192! 193! Sanity check 194! 195 if(n_occ_inact .lt. 0) then 196 call quit('n_occ_inact less than 0') 197 end if 198! 199 if(n_vir_inact .lt. 0) then 200 call quit('n_vir_inact less than 0') 201 end if 202! 203 if(n_gen_inact .lt. 0) then 204 call quit('n_gen_inact less than 0') 205 end if 206! 207 n_occ_aag = n_occ_act*n_occ_act*n_occ_gen 208 n_vir_aag = n_vir_act*n_vir_act*n_vir_gen 209! 210 n_v_a_3 = n_vir_act**3 211! 212 n_occ_integral = n_occ_aag*n_vir_act 213 n_vir_integral = n_vir_aag*n_occ_act 214! 215 if(print_mlcc3 .ge. 3) then 216! 217 write(lupri,*) 218 write(lupri,*) 'Output from mlcc3_active spaces' 219 write(lupri,*) 'n_basis: ', n_basis 220 write(lupri,*) 'n_orbitals: ', n_orbitals 221 write(lupri,*) 'n_occ: ', n_occ 222 write(lupri,*) 'n_occ_act: ', n_occ_act 223 write(lupri,*) 'n_occ_inact: ', n_occ_inact 224 write(lupri,*) 'n_gen_inact: ', n_gen_inact 225 write(lupri,*) 'n_occ_gen: ', n_occ_gen 226 write(lupri,*) 'n_vir: ', n_vir 227 write(lupri,*) 'n_vir_act: ', n_vir_act 228 write(lupri,*) 'n_vir_gen: ', n_vir_gen 229 write(lupri,*) 230! 231 end if 232! 233 end subroutine center_h_import 234! 235! 236 subroutine cvsexci_import 237! 238 implicit none 239! 240 integer :: allocate_ok, i 241! 242 integer, parameter :: maxcore=5, maxion=1 243! 244 integer nrhfcore, irhfcore, nvirion, ivirion 245 logical lcvsexci, lionizexci, lbothexci, lrmcore, lcvsptexci 246! 247 common /cvsepexci/ nrhfcore(8),irhfcore(maxcore,8), & 248 & lcvsexci,lrmcore, & 249 & nvirion(8), ivirion(maxion,8), & 250 & lionizexci, lbothexci, lcvsptexci 251! 252! rmc or cvs 253! 254 if(lcvsexci .or. lrmcore) then 255! 256 cvs_log = lcvsexci 257 rmc_log = lrmcore 258! 259 n_cvs_rmc = nrhfcore(1) 260 261 allocate(cvs_rmc_index(n_cvs_rmc), stat = allocate_ok) 262 if(allocate_ok /= 0) call quit('Something wrong with allocation in center_h_import') 263! 264 do i = 1,n_cvs_rmc 265 cvs_rmc_index(i) = irhfcore(i,1) 266 end do 267! 268 else 269! 270 cvs_log = .false. 271 rmc_log = .false. 272! 273 n_cvs_rmc = 0 274 cvs_rmc_index => null() 275! 276 end if 277! 278! Ionisation 279! 280 if(lionizexci) then 281! 282 ion_log = lionizexci 283! 284 n_ion = nvirion(1) 285 286 allocate(ion_index(n_ion), stat = allocate_ok) 287 if(allocate_ok /= 0) call quit('Something wrong with allocation in center_h_import') 288! 289 do i = 1,n_ion 290 ion_index = ivirion(i,1) 291 end do 292! 293 else 294! 295 ion_log = .false. 296! 297 n_ion = 0 298 ion_index => null() 299! 300 end if 301! 302! Offset indices in case of active space 303! 304 if(mlcc3_active) then 305 if(cvs_log .or. rmc_log) then 306 do i =1,n_cvs_rmc 307 cvs_rmc_index(i) = cvs_rmc_index(i) - n_occ_inact 308 end do 309 end if 310 if(ion_log) then 311 do i =1,n_cvs_rmc 312 ion_index(i) = ion_index(i) - n_vir_inact 313 end do 314 end if 315 end if 316! 317! Sanity check 318! 319 if(cvs_log .or. rmc_log) then 320 do i =1,n_cvs_rmc 321 if(cvs_rmc_index(i) .lt. 1) then 322! 323 write(lupri,*) 'i', i 324 write(lupri,*) 'index', cvs_rmc_index(i) 325 call quit('CVS/RMC index less than 1') 326! 327 else if(cvs_rmc_index(i) .gt. n_occ_act) then 328! 329 write(lupri,*) 'i', i 330 write(lupri,*) 'index', cvs_rmc_index(i) 331 write(lupri,*) 'n_occ_act', n_occ_act 332 call quit('CVS/RMC index greater than n_occ_act') 333! 334 end if 335 end do 336 end if 337! 338 if(ion_log) then 339 do i =1,n_ion 340 if(ion_index(i) .lt. 1) then 341! 342 write(lupri,*) 'i', i 343 write(lupri,*) 'index', cvs_rmc_index(i) 344 call quit('Ion index less than 1') 345! 346 else if(ion_index(i) .gt. n_vir_act) then 347! 348 write(lupri,*) 'i', i 349 write(lupri,*) 'index', ion_index(i) 350 write(lupri,*) 'n_vir_act', n_vir_act 351 call quit('Ion index greater than n_vir_act') 352! 353 end if 354 end do 355 end if 356! 357! 358 end subroutine cvsexci_import 359! 360end module mlcc3_active_spaces 361