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 CCRLXXIM(XIM,ISYXIM,LABEL,LORX,LPDBS,FREQ, 21 & CMO,WORK,LWORK) 22*---------------------------------------------------------------------* 23* 24* Purpose: Calculate the X intermediate needed to calculate 25* the orbital relaxation contributions to second- 26* and higher-order derivatives 27* 28* XIM -- the X intermediate (output) 29* ISYXIM -- symmetry of the X intermediate (input) 30* LABEL -- operator label for the perturbation (input) 31* LORX -- orbital relaxation flag (input) 32* LPDBS -- flag for reorthogonalization contrib. (input) 33* FREQ -- frequency of the perturbation (input) 34* CMO -- HF orbital coefficients 35* 36* Christof Haettig 5-2-1999 37* 38* N.B.: this routine is not yet adapted for QMATP diff. from QMATH 39* 40*---------------------------------------------------------------------* 41 IMPLICIT NONE 42#include "priunit.h" 43#include "ccorb.h" 44#include "ccsdsym.h" 45#include "ccexpfck.h" 46#include "cc1dxfck.h" 47#include "ccfro.h" 48#include "ccroper.h" 49 50 LOGICAL LOCDBG 51 PARAMETER (LOCDBG = .FALSE.) 52 53 INTEGER ISYM0, LUFCK 54 PARAMETER( ISYM0 = 1 ) 55 CHARACTER LABEL0*(8) 56 PARAMETER( LABEL0 = 'HAM0 ' ) 57 58 CHARACTER*(8) LABEL 59 LOGICAL LORX, LPDBS 60 INTEGER ISYXIM, LWORK 61 62#if defined (SYS_CRAY) 63 REAL FREQ, XIM(*), CMO(*), WORK(LWORK) 64 REAL ONE, TWO, ZERO, FTRACE 65#else 66 DOUBLE PRECISION FREQ, XIM(*), CMO(*), WORK(LWORK) 67 DOUBLE PRECISION ONE, TWO, ZERO, FTRACE 68#endif 69 PARAMETER( ONE = 1.0D0, TWO = 2.0D0, ZERO = 0.0D0 ) 70 71 CHARACTER*(10) MODEL 72 LOGICAL NOKAPPA, new 73 INTEGER IFOCK, IADRF, IOPT, IKAPPA, ISYM, IEXPV, IREAL 74 INTEGER KFOCK1,KFOCK0,KOVERLP,KAPPA,KEND1,KSCR1,LWRK1 75 INTEGER IFCK1, KRMAT, KQMATP, KQMATH, IOPER, INDF, INDE1, INDE2 76 INTEGER ISYMA, ISYMI, ISYALP, ISYBET, ISYAL0, ISYBT0 77 INTEGER KCMOQ, KFCKMO, KLAMDPQ, KLAMDHQ, KT1AM, KSCR0 78 INTEGER NBASA, NBASB, NBSA0, NBSB0, NMOA, ISYM1, ISYM2 79 INTEGER KOFF1, KOFF2, KOFF3, KOFF4, KOFF5, KOFF6, KOFF7 80 INTEGER ICMO(8,8), NCMO(8), ICOUNT 81 INTEGER KSAI, KSIA, IORBI, IORBA, IVIRA 82 83* external functions: 84 INTEGER IEFFFOCK 85 INTEGER IEXPECT 86 INTEGER IR1KAPPA 87 INTEGER IROPER 88 INTEGER I1DXFCK 89 90*---------------------------------------------------------------------* 91* some settings: 92*---------------------------------------------------------------------* 93 INDF = 1 94 INDE1 = 1 95 INDE2 = 2 96 97 DO ISYM = 1, NSYM 98 ICOUNT = 0 99 DO ISYM2 = 1, NSYM 100 ISYM1 = MULD2H(ISYM,ISYM2) 101 ICMO(ISYM1,ISYM2) = ICOUNT 102 ICOUNT = ICOUNT + NBAS(ISYM1)*NORBS(ISYM2) 103 END DO 104 NCMO(ISYM) = ICOUNT 105 END DO 106 107*---------------------------------------------------------------------* 108* allocate memory for some intermediates used localy: 109*---------------------------------------------------------------------* 110 KFOCK1 = 1 111 KOVERLP = KFOCK1 + N2BST(ISYXIM) 112 KEND1 = KOVERLP + N2BST(ISYM0) 113 LWRK1 = LWORK - KEND1 114 115 IF (LWRK1 .LT. NBAST) THEN 116 CALL QUIT('Insufficient work space in CCRLXXIM.') 117 END IF 118 119 ! initialize output vector 120 CALL DZERO(XIM,N2BST(ISYXIM)) 121 122 ! read first-order effective Fock matrix from file 123 IFOCK = IEFFFOCK(LABEL,ISYM,1) 124 IEXPV = IEXPECT(LABEL,ISYM) 125 IADRF = IADRFCK(INDF,IFOCK) 126 127 LUFCK = -1 128 CALL WOPEN2(LUFCK,FILFCKEFF,64,0) 129 CALL GETWA2(LUFCK,FILFCKEFF,WORK(KFOCK1),IADRF,N2BST(ISYXIM)) 130 CALL WCLOSE2(LUFCK,FILFCKEFF,'KEEP') 131 132 IF (LOCDBG) THEN 133 FTRACE = ZERO 134 IF (ISYXIM.EQ.1) THEN 135 DO ISYM = 1, NSYM 136 KOFF1 = KFOCK1 + IAODIS(ISYM,ISYM) 137 DO I = 1, NBAS(ISYM) 138 FTRACE = FTRACE + WORK(KOFF1+NBAS(ISYM)*(I-1)+I-1) 139 END DO 140 END DO 141 END IF 142 WRITE (LUPRI,*) 'LABEL:',LABEL 143 WRITE (LUPRI,*) 'ISYXIM,IFOCK,IEXPV:',ISYXIM,IFOCK,IEXPV 144 WRITE (LUPRI,*) 'FTRACE of read matrix:',FTRACE 145 WRITE (LUPRI,*) 'one-electron expect:',EXPVALUE(INDE1,IEXPV) 146 WRITE (LUPRI,*) 'two-electron expect:',EXPVALUE(INDE2,IEXPV) 147 END IF 148 149 ! read zero-order overlap matrix from file and square up 150 CALL RDONEL('OVERLAP ',.TRUE.,WORK(KEND1),NBAST) 151 CALL CCSD_SYMSQ(WORK(KEND1),ISYM0,WORK(KOVERLP)) 152 153 154 ! transform effective Fock matrix to MO basis 155 CALL CC_EFFCKMO(WORK(KFOCK1),ISYXIM,CMO,WORK(KOVERLP), 156 & WORK(KEND1),LWRK1) 157 158 IF (LOCDBG) THEN 159 FTRACE = ZERO 160 IF (ISYXIM.EQ.1) THEN 161 DO ISYM = 1, NSYM 162 KOFF1 = KFOCK1 + IAODIS(ISYM,ISYM) 163 DO I = 1, NBAS(ISYM) 164 FTRACE = FTRACE + WORK(KOFF1+NBAS(ISYM)*(I-1)+I-1) 165 END DO 166 END DO 167 END IF 168 WRITE (LUPRI,*) 'LABEL:',LABEL 169 WRITE (LUPRI,*) 'ISYXIM:',ISYXIM 170 WRITE (LUPRI,*) 'FTRACE of matrix generated in CC_EFCKMO:', 171 & FTRACE 172 END IF 173 174 ! add 2 times to X intermediate 175 CALL DAXPY(N2BST(ISYXIM),TWO,WORK(KFOCK1),1,XIM,1) 176 177 IF (LOCDBG) THEN 178 WRITE (LUPRI,*) 179 & 'CCRLXXIM> direct contribution to X intermediate:' 180 CALL CC_PRONELAO(XIM,ISYXIM) 181 END IF 182 183 184 IF (LORX .OR. LPDBS) THEN 185 186 IOPER = IROPER(LABEL,ISYM) 187 IREAL = ISYMAT(IOPER) 188 IFCK1 = I1DXFCK('HAM0 ','R1 ',LABEL,FREQ,ISYM) 189 190 KFOCK0 = KEND1 191 KFOCK1 = KFOCK0 + N2BST(ISYM0) 192 KAPPA = KFOCK1 + N2BST(ISYXIM) 193 KRMAT = KAPPA + 2*NALLAI(ISYXIM) 194 KQMATP = KRMAT + N2BST(ISYXIM) 195 KQMATH = KQMATP + N2BST(ISYXIM) 196 KSCR1 = KQMATH + N2BST(ISYXIM) 197 KEND1 = KSCR1 + N2BST(ISYXIM) 198 LWRK1 = LWORK - KEND1 199 200 IF (LWRK1 .LT. 0) THEN 201 CALL QUIT('Insufficient work space in CCRLXXIM.') 202 END IF 203 204 IF (LORX) THEN 205 IKAPPA = IR1KAPPA(LABEL,FREQ,ISYM) 206 CALL CC_RDHFRSP('R1 ',IKAPPA,ISYM,WORK(KAPPA)) 207 ELSE 208 CALL DZERO(WORK(KAPPA),2*NALLAI(ISYXIM)) 209 END IF 210 211 CALL CC_GET_RMAT(WORK(KRMAT),IOPER,1,ISYXIM,WORK(KEND1),LWRK1) 212 213 new = .true. 214 if (.not. new) then 215 CALL DSCAL(N2BST(ISYXIM),-1.0d0,WORK(KRMAT),1) 216 end if 217 218 NOKAPPA = .FALSE. 219 CALL CC_QMAT(WORK(KQMATP),WORK(KQMATH),WORK(KRMAT),WORK(KAPPA), 220 & IREAL,ISYXIM,NOKAPPA,CMO,WORK(KEND1),LWRK1) 221 222 if (new) then 223 DO ISYM1 = 1, NSYM 224 ISYM2 = MULD2H(ISYM1,ISYXIM) 225 KOFF1 = KQMATP + IAODIS(ISYM1,ISYM2) 226 KOFF2 = KQMATH + IAODIS(ISYM2,ISYM1) 227 CALL TRSREC(NBAS(ISYM1),NBAS(ISYM2), 228 & WORK(KOFF1),WORK(KOFF2)) 229 END DO 230 CALL DSCAL(N2BST(ISYXIM),-ONE,WORK(KQMATH),1) 231 end if 232 233 IFOCK = IEFFFOCK(LABEL0,ISYM,1) 234 IEXPV = IEXPECT(LABEL0,ISYM) 235 IADRF = IADRFCK(INDF,IFOCK) 236 CALL WOPEN2(LUFCK,FILFCKEFF,64,0) 237 CALL GETWA2(LUFCK,FILFCKEFF,WORK(KFOCK0),IADRF,N2BST(ISYM0)) 238 CALL WCLOSE2(LUFCK,FILFCKEFF,'KEEP') 239 240 IF (LOCDBG) THEN 241 FTRACE = ZERO 242 DO ISYM = 1, NSYM 243 KOFF1 = KFOCK0 + IAODIS(ISYM,ISYM) 244 DO I = 1, NBAS(ISYM) 245 FTRACE = FTRACE + WORK(KOFF1+NBAS(ISYM)*(I-1)+I-1) 246 END DO 247 END DO 248 WRITE (LUPRI,*) 'LABEL:',LABEL0 249 WRITE (LUPRI,*) 'ISYXIM,IFOCK,IEXPV:',ISYXIM,IFOCK,IEXPV 250 WRITE (LUPRI,*) 'FTRACE of read matrix :',FTRACE 251 WRITE (LUPRI,*) 'one-electron expect:',EXPVALUE(INDE1,IEXPV) 252 WRITE (LUPRI,*) 'two-electron expect:',EXPVALUE(INDE2,IEXPV) 253 END IF 254 255 CALL CC_EFFCKMO(WORK(KFOCK0),ISYM0,CMO,WORK(KOVERLP), 256 & WORK(KEND1),LWRK1) 257 258 IF (LOCDBG) THEN 259 FTRACE = ZERO 260 DO ISYM = 1, NSYM 261 KOFF1 = KFOCK0 + IAODIS(ISYM,ISYM) 262 DO I = 1, NBAS(ISYM) 263 FTRACE = FTRACE + WORK(KOFF1+NBAS(ISYM)*(I-1)+I-1) 264 END DO 265 END DO 266 WRITE (LUPRI,*) 'LABEL:',LABEL0 267 WRITE (LUPRI,*) 'ISYXIM:',ISYM0 268 WRITE (LUPRI,*) 269 & 'FTRACE of matrix generated by CC_EFFCKMO:',FTRACE 270 END IF 271 272 ! calculate F^0 x Q, result is put into scr 273 CALL CC_MMOMMO('N','N',ONE,WORK(KFOCK0),ISYM0, 274 & WORK(KQMATH),ISYXIM,ZERO,WORK(KSCR1),ISYXIM) 275 276 IF (LOCDBG) THEN 277 WRITE (LUPRI,*) 278 & 'CCRLXXIM> relax. contrib. 1 to X intermediate:' 279 CALL CC_PRONELAO(WORK(KSCR1),ISYXIM) 280 FTRACE = ZERO 281 DO ISYM = 1, NSYM 282 KOFF1 = KSCR1 + IAODIS(ISYM,ISYM) 283 DO I = 1, NBAS(ISYM) 284 FTRACE = FTRACE + WORK(KOFF1+NBAS(ISYM)*(I-1)+I-1) 285 END DO 286 END DO 287 WRITE (LUPRI,*) 'trace:',FTRACE 288 END IF 289 290 ! add to result matrix: 291 CALL DAXPY(N2BST(ISYXIM),ONE,WORK(KSCR1),1,XIM,1) 292 293 ! read contribution from 'one-index' transformed density 294 ! from file 295 IADRF = IADR1DXF(1,IFCK1) 296 CALL WOPEN2(LUFCK,FIL1DXFCK,64,0) 297 CALL GETWA2(LUFCK,FIL1DXFCK,WORK(KFOCK1),IADRF,N2BST(ISYXIM)) 298 CALL WCLOSE2(LUFCK,FIL1DXFCK,'KEEP') 299 300 CALL CC_EFFCKMO(WORK(KFOCK1),ISYXIM,CMO,WORK(KOVERLP), 301 & WORK(KEND1),LWRK1) 302 303 IF (LOCDBG) THEN 304 WRITE (LUPRI,*) 305 & 'CCRLXXIM> relax. contrib. 2 to X intermediate:' 306 CALL CC_PRONELAO(WORK(KFOCK1),ISYXIM) 307 FTRACE = ZERO 308 DO ISYM = 1, NSYM 309 KOFF1 = KFOCK1 + IAODIS(ISYM,ISYM) 310 DO I = 1, NBAS(ISYM) 311 FTRACE = FTRACE + WORK(KOFF1+NBAS(ISYM)*(I-1)+I-1) 312 END DO 313 END DO 314 WRITE (LUPRI,*) 'trace:',FTRACE 315 END IF 316 317 ! add to result matrix: 318 CALL DAXPY(N2BST(ISYXIM),ONE,WORK(KFOCK1),1,XIM,1) 319 END IF 320 321 IF (LOCDBG) THEN 322 WRITE (LUPRI,*) 'CCRLXXIM> final result for X intermediate:' 323 CALL CC_PRONELAO(XIM,ISYXIM) 324 FTRACE = ZERO 325 DO ISYM = 1, NSYM 326 KOFF1 = IAODIS(ISYM,ISYM) 327 DO I = 1, NBAS(ISYM) 328 FTRACE = FTRACE + XIM(KOFF1+NBAS(ISYM)*(I-1)+I) 329 END DO 330 END DO 331 WRITE (LUPRI,*) 'trace of X intermediate:',FTRACE 332 END IF 333 334 RETURN 335 END 336*======================================================================* 337