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! 18C 19*=====================================================================* 20 SUBROUTINE CCREO2CON(IREO,NREO,IRDOTS,RCONS, 21 & MXVEC,WORK,LWORK) 22*---------------------------------------------------------------------* 23* 24* Purpose: calculate contributions from second-order orbital 25* reorthogonalization and the coupling of reorthog. 26* and relaxation to response functions involving 27* perturbation-dependent basis sets. 28* 29* IREO -- array with the perturbation indeces 30* NREO -- length of IREO 31* IRDOTS -- matrix containing the second indeces 32* RCONS -- matrix with the contributions 33* MXVEC -- leading dimension of IRDOTS and RCONS 34* 35* Christof Haettig 11-6-1999 36* 37*---------------------------------------------------------------------* 38 IMPLICIT NONE 39#include "priunit.h" 40#include "dummy.h" 41#include "ccorb.h" 42#include "ccsdsym.h" 43#include "ccexpfck.h" 44#include "cc1dxfck.h" 45#include "ccr1rsp.h" 46#include "ccfro.h" 47#include "ccroper.h" 48#include "iratdef.h" 49#include "inftap.h" 50 51 LOGICAL LOCDBG 52 PARAMETER (LOCDBG = .FALSE.) 53 54 INTEGER ISYM0, LUFCK 55 PARAMETER( ISYM0 = 1) 56 CHARACTER LABEL0*(8) 57 PARAMETER( LABEL0 = 'HAM0 ' ) 58 59 LOGICAL LORX, LFOCK0 60 INTEGER LWORK, NREO, MXVEC 61 62 INTEGER IREO(NREO), IRDOTS(MXVEC,NREO) 63 64#if defined (SYS_CRAY) 65 REAL TWO, FREQ, RCONS(MXVEC,NREO), WORK(LWORK) 66#else 67 DOUBLE PRECISION TWO, FREQ, RCONS(MXVEC,NREO), WORK(LWORK) 68#endif 69 PARAMETER (TWO = 2.0D0) 70 71 CHARACTER*(10) MODEL 72 CHARACTER*(8) LABELA, LABELB 73 LOGICAL LORXA, LORXB, LPDBSA, LPDBSB, NOKAPPA 74 INTEGER NCMOT, NASHT, N2ASHX, LCINDX, IDXREO, IKAPPA 75 INTEGER KCMO, KFOCK0, KOVERLP, KEND1, LWRK1, IFOCK, IADRF, ISYM 76 INTEGER ISYMA, ISYMB, IOPERA, IOPERB, KR2EFF, KRMATA, KAPASQ 77 INTEGER KAPPAA, KAPBSQ, KAPPAB, KSCR1, KEND2, LWRK2, KRMATB 78 INTEGER KQMATPA, KQMATPB 79 INTEGER IOPT, IVEC, IKAPPB, ISAMA, ISAMB 80 81* external functions: 82 INTEGER IEFFFOCK 83 INTEGER ILSTSYM, ILSTSYMRLX 84 INTEGER IROPER 85#if defined (SYS_CRAY) 86 REAL DDOT 87#else 88 DOUBLE PRECISION DDOT 89#endif 90 91*---------------------------------------------------------------------* 92* check, if there is anything at all to do: 93*---------------------------------------------------------------------* 94 95 IF (NREO.LE.0) RETURN 96 97*---------------------------------------------------------------------* 98* get some constants from sirius common block: 99*---------------------------------------------------------------------* 100 101 CALL CC_SIRINF(NCMOT,NASHT,N2ASHX,LCINDX) 102 103*---------------------------------------------------------------------* 104* allocate memory for perturbation-independent stuff needed: 105*---------------------------------------------------------------------* 106 KCMO = 1 107 KFOCK0 = KCMO + NCMOT 108 KOVERLP = KFOCK0 + N2BST(ISYM0) 109 KEND1 = KOVERLP + N2BST(ISYM0) 110 LWRK1 = LWORK - KEND1 111 112 IF (LWRK1.LT.0) THEN 113 CALL QUIT('Insufficient work space in CCREO2CON.') 114 END IF 115 116*---------------------------------------------------------------------* 117* read MO coefficients from file: 118*---------------------------------------------------------------------* 119 120 CALL GPOPEN(LUSIFC,'SIRIFC','OLD',' ','UNFORMATTED',IDUMMY, 121 & .FALSE.) 122 REWIND LUSIFC 123 CALL MOLLAB('SIR IPH ',LUSIFC,LUPRI) 124 READ (LUSIFC) 125 READ (LUSIFC) 126 CALL READI(LUSIFC,IRAT*NCMOT,WORK(KCMO)) 127 CALL GPCLOSE(LUSIFC,'KEEP') 128 129*---------------------------------------------------------------------* 130* read overlap matrix from file: 131*---------------------------------------------------------------------* 132 133 IF (LWRK1.LT.NBAST) THEN 134 CALL QUIT('Insufficient work space in CCREO2CON.') 135 END IF 136 137 CALL RDONEL('OVERLAP ',.TRUE.,WORK(KEND1),NBAST) 138 CALL CCSD_SYMSQ(WORK(KEND1),ISYM0,WORK(KOVERLP)) 139 140*---------------------------------------------------------------------* 141* read zeroth-order effective Fock matrix if available: 142*---------------------------------------------------------------------* 143 IFOCK = IEFFFOCK(LABEL0,ISYM,1) 144 145 IF (LEXPFCK(2,IFOCK)) THEN 146 IADRF = IADRFCK(1,IFOCK) 147 148 LUFCK = -1 149 CALL WOPEN2(LUFCK,FILFCKEFF,64,0) 150 CALL GETWA2(LUFCK,FILFCKEFF,WORK(KFOCK0),IADRF,N2BST(ISYM0)) 151 CALL WCLOSE2(LUFCK,FILFCKEFF,'KEEP') 152 153 CALL CC_EFFCKMO(WORK(KFOCK0),ISYM0,WORK(KCMO),WORK(KOVERLP), 154 & WORK(KEND1),LWRK1) 155 156 LFOCK0 = .TRUE. 157 ELSE 158 LFOCK0 = .FALSE. 159 END IF 160 161*---------------------------------------------------------------------* 162* start loop over first index: 163*---------------------------------------------------------------------* 164 DO IDXREO = 1, NREO 165 166 IKAPPA = IREO(IDXREO) 167 ISYMA = ILSTSYMRLX('R1',IKAPPA) 168 LABELA = LRTHFLBL(IKAPPA) 169 IOPERA = IROPER(LABELA,ISYMA) 170 ISAMA = ISYMAT(IOPERA) 171 LORXA = .TRUE. 172 LPDBSA = LPDBSOP(IOPERA) 173 174 IF (LOCDBG) THEN 175 WRITE (LUPRI,*) 'CCREO2CON> IDXREO:',IDXREO 176 WRITE (LUPRI,*) 'CCREO2CON> IKAPPA:',IKAPPA 177 WRITE (LUPRI,*) 'CCREO2CON> LABELA:',LABELA 178 WRITE (LUPRI,*) 'CCREO2CON> LORXA :',LORXA 179 WRITE (LUPRI,*) 'CCREO2CON> LPDBSA:',LPDBSA 180 END IF 181 182 KR2EFF = KEND1 183 KQMATPA = KR2EFF + N2BST(ISYM0) 184 KRMATA = KQMATPA + N2BST(ISYMA) 185 KAPASQ = KRMATA + N2BST(ISYMA) 186 KAPPAA = KAPASQ + N2BST(ISYMA) 187 KQMATPB = KAPPAA + 2*NALLAI(ISYMA) 188 KRMATB = KQMATPB + N2BST(ISYMA) 189 KAPBSQ = KRMATB + N2BST(ISYMA) 190 KAPPAB = KAPBSQ + N2BST(ISYMA) 191 KSCR1 = KAPPAB + NALLAI(ISYMA) 192 KEND2 = KSCR1 + N2BST(ISYMA) 193 194 LWRK2 = LWORK - KEND2 195 196 IF (LWRK2 .LT. 0) THEN 197 CALL QUIT('Insufficient memory in CCREO2CON.') 198 END IF 199 200 IF (LORXA) THEN 201 CALL CC_RDHFRSP('R1 ',IKAPPA,ISYMA,WORK(KAPPAA)) 202 CALL CCKAPPASQ(WORK(KAPASQ),WORK(KAPPAA),ISYMA,'N') 203 END IF 204 205 IF (LPDBSA) THEN 206 CALL CC_GET_RMAT(WORK(KSCR1),IOPERA,1,ISYMA, 207 & WORK(KEND2),LWRK2) 208 209 NOKAPPA = .TRUE. 210 CALL CC_QMAT(WORK(KQMATPA),WORK(KRMATA),WORK(KSCR1),DUMMY, 211 & ISAMA,ISYMA,NOKAPPA,WORK(KCMO),WORK(KEND2),LWRK2) 212 END IF 213 214*---------------------------------------------------------------------* 215* loop over second index: 216*---------------------------------------------------------------------* 217 IVEC = 1 218 219 DO WHILE (IRDOTS(IVEC,IDXREO).NE.0 .AND. IVEC.LE.MXVEC) 220 221C WRITE (LUPRI,*) 'CCREO2CON> IVEC = ',IVEC 222 223 IKAPPB = IRDOTS(IVEC,IDXREO) 224 ISYMB = ILSTSYMRLX('R1',IKAPPB) 225 LABELB = LRTHFLBL(IKAPPB) 226 IOPERB = IROPER(LABELB,ISYMB) 227 ISAMB = ISYMAT(IOPERB) 228 LORXB = .TRUE. 229 LPDBSB = LPDBSOP(IOPERB) 230 231 IF (LOCDBG) THEN 232 WRITE (LUPRI,*) 'CCREO2CON> IVEC :',IVEC 233 WRITE (LUPRI,*) 'CCREO2CON> IKAPPB:',IKAPPB 234 WRITE (LUPRI,*) 'CCREO2CON> LABELB:',LABELB 235 WRITE (LUPRI,*) 'CCREO2CON> LORXB :',LORXB 236 WRITE (LUPRI,*) 'CCREO2CON> LPDBSB:',LPDBSB 237 END IF 238 239 IF (ISYMB.NE.ISYMA) THEN 240 WRITE (LUPRI,*) IDXREO, ISYMA, IVEC, ISYMB 241 CALL QUIT('symmetry mismatch in CCREO2CON.') 242 END IF 243 244 IF (LORXB) THEN 245 CALL CC_RDHFRSP('R1 ',IKAPPB,ISYMB,WORK(KAPPAB)) 246 CALL CCKAPPASQ(WORK(KAPBSQ),WORK(KAPPAB),ISYMB,'N') 247 END IF 248 249 IF (LPDBSB) THEN 250 CALL CC_GET_RMAT(WORK(KSCR1),IOPERB,1,ISYMB, 251 & WORK(KEND2),LWRK2) 252 253 NOKAPPA = .TRUE. 254 CALL CC_QMAT(WORK(KQMATPB),WORK(KRMATB),WORK(KSCR1),DUMMY, 255 & ISAMB,ISYMB,NOKAPPA,WORK(KCMO),WORK(KEND2),LWRK2) 256 END IF 257 258 CALL DZERO(WORK(KR2EFF),N2BST(ISYM0)) 259 260 IF (LORXB.AND.LPDBSA) THEN 261 CALL CC_MMOMMO('N','N',+1.0D0,WORK(KAPBSQ),ISYMB, 262 & WORK(KRMATA),ISYMA,1.0D0,WORK(KR2EFF),1) 263 CALL CC_MMOMMO('N','N',-1.0D0,WORK(KRMATA),ISYMA, 264 & WORK(KAPBSQ),ISYMB,1.0D0,WORK(KR2EFF),1) 265 END IF 266 267 IF (LORXA.AND.LPDBSB) THEN 268 CALL CC_MMOMMO('N','N',+1.0D0,WORK(KAPASQ),ISYMA, 269 & WORK(KRMATB),ISYMB,1.0D0,WORK(KR2EFF),1) 270 CALL CC_MMOMMO('N','N',-1.0D0,WORK(KRMATB),ISYMB, 271 & WORK(KAPASQ),ISYMA,1.0D0,WORK(KR2EFF),1) 272 END IF 273 274 IF (LPDBSA.AND.LPDBSB) THEN 275 CALL QUIT('CCREO2CON NOT YET IMPLEMENTED '// 276 & 'FOR 2. DERIVATIVES.') 277 END IF 278 279 RCONS(IVEC,IDXREO) = TWO*DDOT(N2BST(ISYM0),WORK(KR2EFF),1, 280 & WORK(KFOCK0),1) 281C WRITE (LUPRI,*) 'CCREO2CON>',RCONS(IVEC,IDXREO) 282 283 IVEC = IVEC + 1 284 END DO 285 286 END DO 287 288*---------------------------------------------------------------------* 289* print the results: 290*---------------------------------------------------------------------* 291 IF (LOCDBG) THEN 292 WRITE(LUPRI,*) 293 & 'CCREO2CON> results for X intermediate contribs.:' 294 IF (MXVEC.NE.0) THEN 295 DO IDXREO = 1, NREO 296 WRITE (LUPRI,*) 'IDXREO = ',IDXREO 297 IVEC = 1 298 DO WHILE(IRDOTS(IVEC,IDXREO).NE.0 .AND. IVEC.LE.MXVEC) 299 WRITE(LUPRI,'(A,2I5,2X,E19.12)') 'CCREO2CON> ', 300 & IREO(IDXREO),IRDOTS(IVEC,IDXREO),RCONS(IVEC,IDXREO) 301 IVEC = IVEC + 1 302 END DO 303 END DO 304 ELSE 305 WRITE (LUPRI,*) 'MXVEC.EQ.0 --> nothing calculated.' 306 END IF 307 END IF 308 309 RETURN 310 END 311*======================================================================* 312