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! look for HJMAERKE; (1) symort of UMAT should be done 21! 22!961213-hjaaj: 23!AVESET: disable SUPSYM if no supersymmetry found 24! HJMAERKE: this means that debug SUPSYM for abelian group 25! cannot be performed without commenting new code in AVESET! 26!940815-hjaaj: 27!AVESE1,AVEDEG,AVEGTU,AVEPHA: use THRSSY from INFINP (new in INFINP) 28! instead of THREQL parameter; include INFINP 29!AVESE1: Check and reset THRSSY if needed; more print if disentangling 30!940427-hjaaj: SCALAR directives for Cray (as SCALSUB does not work 31! for Cray because directive has to be after SUBROUTINE statement) 32!930715-hjaaj 33!AVESE0: define MXDGSS and NINFSS(*,3) 34!AVEDEG: define MXDGSS 35!=========================================================================== 36C /* Deck averag */ 37 SUBROUTINE AVERAG (VECS,NVDIM,NVSIM) 38C 39C 7-Jul-1992 hjaaj+hh 40C 41#include "implicit.h" 42#include "priunit.h" 43 DIMENSION VECS(NVDIM,NVSIM) 44C 45C Used from common blocks: 46C INFINP : SUPSYM 47C INFVAR : NWOPT,JWOPSY 48C 49#include "maxorb.h" 50#include "infinp.h" 51#include "infvar.h" 52#include "infpri.h" 53C 54 CALL QENTER('AVERAG') 55 IF (.NOT. SUPSYM) GO TO 9999 56 IF (NWOPT .EQ. 0) GO TO 9999 57 IF (JWOPSY .NE. 1) GO TO 910 58 IF (NVDIM .LT. NWOPT) GO TO 920 59 IF (SUPSYM) THEN 60 CALL AVERSS(VECS,NVDIM,NVSIM) 61 END IF 62 9999 CALL QEXIT('AVERAG') 63 RETURN 64C 65C *** Error sections 66C 67 910 CONTINUE 68 WRITE (LUERR,9010) JWOPSY 69 CALL QTRACE(LUERR) 70 CALL QUIT('AVERAG ERROR, operator symmetry JWOPSY .ne. 1') 71 9010 FORMAT(/' ERROR-AVERAG, operator symmetry .ne. 1; JWOPSY =',I8) 72C 73 920 CONTINUE 74 WRITE (LUERR,9020) NVDIM,NWOPT 75 CALL QTRACE(LUERR) 76 CALL QUIT('AVERAG ERROR, vector length lt NWOPT') 77 9020 FORMAT(/' ERROR-AVERAG, vector length',I8,' is less than NWOPT', 78 * I8) 79C 80 END 81C /* Deck averss */ 82 SUBROUTINE AVERSS (VECS,NVDIM,NVSIM) 83C 84C 7-Jul-1992 Hans Joergen Aa. Jensen + Hinne Hettema 85C 86C Purpose: Average degenerate super symmetries 87C (which are different components of a degenerate irrep). 88C This corresponds to doing a microcanonical average. 89C 90#include "implicit.h" 91 DIMENSION VECS(NVDIM,NVSIM) 92C 93C Used from common blocks: 94C INFORB : NSSYM,NORBSS(),IORBSS(),NINFSS(*,3), NOCCT,NORBT 95C INFIND : ISW(*),ISSORD(*) 96C INFVAR : NWOPT 97C 98#include "maxash.h" 99#include "maxorb.h" 100#include "priunit.h" 101#include "infind.h" 102#include "inforb.h" 103#include "infvar.h" 104#include "infpri.h" 105C 106 PARAMETER (MAXDEG = 20, D1 = 1.0D0) 107 DIMENSION ISSDEG(MAXDEG),ISSOFF(MAXDEG),IDEGV(MAXDEG) 108C 109Chjaaj-may2000: 110C Dynamic allocation of KLWOP with adjustable run-time dimensions 111 DIMENSION KLWOP(NOCCT,NORBT) 112C 113 CALL QENTER('AVERSS') 114 CALL MAKE_KLWOP(KLWOP) 115 IF (IPRAVE .GT. 100) THEN 116 WRITE(LUPRI,'(/A)') ' AVERSS, Incoming orbital vector(s):' 117 CALL OUTPUT(VECS,1,NWOPT,1,NVSIM,NVDIM,NVSIM,-1,LUPRI) 118 END IF 119C 120 DO 800 ISSYMR = 1,NSSYM 121C Skip this supsym if not degenerate or not "root supsym" 122 IF (NINFSS(ISSYMR,1) .EQ. 1) GO TO 800 123 IF (NINFSS(ISSYMR,2) .NE. ISSYMR) GO TO 800 124 NDEG = 0 125 DO 100 ISSYM = ISSYMR,NSSYM 126 IF (NINFSS(ISSYM,2) .EQ. ISSYMR) THEN 127 NDEG = NDEG + 1 128 ISSDEG(NDEG) = ISSYM 129 ISSOFF(NDEG) = IORBSS(ISSYM) - IORBSS(ISSYMR) 130 END IF 131 100 CONTINUE 132 IF (IPRAVE .GT. 5 .OR. NDEG .NE. NINFSS(ISSYMR,1)) THEN 133 WRITE(LUPRI,'(/A,I3,A/A,(T35,8I5))') 134 & ' AVERSS averaging ',NDEG,' degenerate super symmetries', 135 & ' namely sup.sym.s no.',(ISSDEG(IDEG),IDEG=1,NDEG) 136 END IF 137 IF (NDEG .NE. NINFSS(ISSYMR,1)) THEN 138 WRITE(LUPRI,'(/A)') 139 & ' AVERSS-software error: inconsistent degeneracy in NINFSS' 140 CALL AVEDMP(LUPRI) 141 CALL QUIT('AVERSS: inconsistent degeneracy in NINFSS') 142 END IF 143C 144 AVFAC = NDEG 145 AVFAC = D1 / AVFAC 146C 147 IOFFSR = IORBSS(ISSYMR) 148 DO 480 IMOSS = IOFFSR+1,IOFFSR+NORBSS(ISSYMR) 149 IMO = ISSORD(IMOSS) 150 IMOW = ISW(IMO) 151 IF (IMOW .GT. NOCCT) GO TO 480 152 DO 470 JMOSS = IMOSS+1,IOFFSR+NORBSS(ISSYMR) 153 JMO = ISSORD(JMOSS) 154 JMOW = ISW(JMO) - NISHT 155 IF (JMOW .LE. 0) GO TO 470 156Chjaaj ... this is an inactive orbital 157 IF (KLWOP(IMOW,JMOW) .EQ. 0) GO TO 470 158 IDEGV(1) = KLWOP(IMOW,JMOW) 159 NERR = 0 160 DO 410 IDEG = 2,NDEG 161 KMO = ISSORD(IMOSS + ISSOFF(IDEG)) 162 LMO = ISSORD(JMOSS + ISSOFF(IDEG)) 163 LMOW = ISW(LMO) - NISHT 164 IF (LMOW .LE. 0) THEN 165 IDEGV(IDEG) = 0 166 NERR = NERR + 1 167 WRITE(LUPRI,'(/A)') 168 & 'AVERSS error: degenerate orbitals not same orbital class.' 169 ELSE 170 IDEGV(IDEG) = KLWOP( ISW(KMO) , LMOW ) 171 IF (IDEGV(IDEG) .EQ. 0) NERR = NERR + 1 172 END IF 173 410 CONTINUE 174 IF (IPRAVE .GT. 10 .OR. NERR .GT. 0) THEN 175 IF (NERR .GT. 0) THEN 176 WRITE(LUPRI,'(/A/A/)') 177 & ' AVERSS error: not all degenerate rotations included', 178 & ' Probable error: some but not all'// 179 & ' of the degenerate orbitals have been frozen.' 180 END IF 181 WRITE(LUPRI,'(A,(T35,4I10))') 182 & ' AVERSS averaging orb.rot. no.s', 183 & (IDEGV(IDEG),IDEG=1,NDEG) 184 IF (NERR .GT. 0) CALL AVEDMP(LUPRI) 185 END IF 186 IF (NERR .GT. 0) THEN 187 CALL QTRACE(LUPRI) 188 CALL QUIT( 189 & 'AVERSS error: not all degenerate rotations included') 190 END IF 191 DO 460 IVSIM = 1,NVSIM 192 XKAPIJ = VECS(IDEGV(1),IVSIM) 193 DO 420 IDEG = 2,NDEG 194 XKAPIJ = XKAPIJ + VECS(IDEGV(IDEG),IVSIM) 195 420 CONTINUE 196 XKAPIJ = AVFAC*XKAPIJ 197 DO 430 IDEG = 1,NDEG 198 VECS(IDEGV(IDEG),IVSIM) = XKAPIJ 199 430 CONTINUE 200 460 CONTINUE 201 470 CONTINUE 202 480 CONTINUE 203 800 CONTINUE 204C 205 IF (IPRAVE .GT. 100) THEN 206 WRITE(LUPRI,'(/A)') ' AVERSS, Averaged orbital vector(s):' 207 CALL OUTPUT(VECS,1,NWOPT,1,NVSIM,NVDIM,NVSIM,-1,LUPRI) 208 END IF 209C 210 CALL QEXIT('AVERSS') 211 RETURN 212 END 213C /* Deck aveset */ 214 SUBROUTINE AVESET(CMO,WRK,KFREE,LFREE) 215C 216C 2-Jul-1992 hjaaj : core allocation for AVESE1 217C 218#include "implicit.h" 219 REAL*8 CMO(*), WRK(*) 220C 221C Used from common blocks: 222C INFINP: SUPSYM 223C INFORB: N2ORBX, MXSSYM 224C 225#include "maxorb.h" 226#include "priunit.h" 227#include "infinp.h" 228#include "inforb.h" 229#include "infpri.h" 230C 231 CALL QENTER('AVESET') 232 IF (SUPSYM) THEN 233 IF (IPRAVE .GE. 3) THEN 234 WRITE(LUPRI,'(//A)') ' TEST OUTPUT FROM AVESET' 235 END IF 236 KFRSAV = KFREE 237 LFRSAV = LFREE 238 CALL MEMGET('REAL',KTKMO,N2ORBX,WRK,KFREE,LFREE) 239 CALL MEMGET('INTE',KNDEGS,MXSSYM,WRK,KFREE,LFREE) 240 CALL AVESE1(CMO,WRK(KTKMO),WRK(KNDEGS), 241 * WRK,KFREE,LFREE) 242 CALL MEMREL('AVESET',WRK,KFRSAV,KFRSAV,KFREE,LFREE) 243C 244C Check if supersymmetry is identical to "D2h" symmetry 245C 246 IF (NSSYM .LE. NSYM .AND. MXDGSS .EQ. 1) THEN 247 DO 200 ISYM = 1,NSYM 248 JSYM = 0 249 DO 100 ISSYM = 1,NSSYM 250 IF (NINFSS(ISSYM,3) .EQ. ISYM) JSYM = JSYM + 1 251 100 CONTINUE 252 IF (JSYM .GT. 1) GO TO 300 253 200 CONTINUE 254 IF (IPRAVE .GE. 1) THEN 255 WRITE(LUPRI,'(/A)')' AVESET: No "super symmetry" found.' 256 END IF 257 SUPSYM = .FALSE. 258 GO TO 300 259 END IF 260 300 CONTINUE 261 END IF 262 IF (.NOT.SUPSYM) THEN 263 CALL AVESE0 264 END IF 265 CALL QEXIT('AVESET') 266 RETURN 267 END 268C /* Deck avese0 */ 269 SUBROUTINE AVESE0 270C 271C 7-Jul-1992 hjaaj 272C Set super symmetry information as symmetry information 273C (when super symmetry is not to be used() 274C 275#include "implicit.h" 276C 277C Used from common blocks: 278C INFORB : NSSYM, NORBSS(), IORBSS(); NSYM, NORB(), IORB(), NORBT 279C INFIND : ISSMO(), ISSORD(); ISMO() 280C 281#include "maxash.h" 282#include "maxorb.h" 283#include "inforb.h" 284#include "infind.h" 285C 286 CALL QENTER('AVESE0') 287 NSSYM = NSYM 288 DO 100 ISYM = 1,NSYM 289 NORBSS(ISYM) = NORB(ISYM) 290 IORBSS(ISYM) = IORB(ISYM) 291 NINFSS(ISYM,1) = 1 292 NINFSS(ISYM,2) = ISYM 293 NINFSS(ISYM,3) = ISYM 294 100 CONTINUE 295 DO 200 IMO = 1,NORBT 296 ISSMO(IMO) = ISMO(IMO) 297 ISSORD(IMO) = IMO 298 200 CONTINUE 299 MXDGSS = 1 300 CALL QEXIT('AVESE0') 301 RETURN 302 END 303C /* Deck avese1 */ 304 SUBROUTINE AVESE1(CMO,TKMO,NDEGSS,WRK,KFREE,LFREE) 305C 306C Written 29-Jun-1992 Hans Joergen Aa. Jensen and Hinne Hettema 307C 308C Purpose: 309C Find super symmetry and degeneracy of molecular orbitals in CMO 310C Disentangle mixed degenerate orbitals 311C Coordinate relative phases of degenerate orbitals 312C 313C Input: 314C CMO; molecular orbitals 315C 316C Output: 317C 318#include "implicit.h" 319 REAL*8 CMO(NCMOT), TKMO(NORBT,NORBT), WRK(*) 320 INTEGER NDEGSS(*) ! DIMENSION NDEGSS(MXSSYM) 321C 322#include "threql.h" 323 PARAMETER (TMXSSY = 1.0D-4) 324C 325C Used from common blocks: 326C INFINP : THRSSY 327C INFORB : NSYM,NSSYM,NORBT,NORB(),? 328C INFIND : ISSMO(),ISSORD(),? 329C INFPRI : IPRAVE,NINFO,? 330C 331#include "maxash.h" 332#include "maxorb.h" 333#include "priunit.h" 334#include "infinp.h" 335#include "inforb.h" 336#include "infind.h" 337#include "infpri.h" 338C 339C 340 CALL QENTER('AVESE1') 341 IF (THRSSY .GT. TMXSSY) THEN 342 NWARN = NWARN + 1 343 WRITE (LUPRI,'(/A/2(A,1P,D10.2))') ' AVESE1-WARNING: SUPSYM'// 344 & ' threshold for degeneracy and "zero" elements too big;', 345 & ' value reset from',THRSSY,' to',TMXSSY 346 THRSSY = TMXSSY 347 END IF 348 IF (THRSSY .LT. THREQL) THEN 349C ... THREQL must be .gt. accuracy of diagonalization routine 350 NINFO = NINFO + 1 351 WRITE (LUPRI,'(/A/2(A,1P,D10.2))') ' AVESE1-INFO:'// 352 & ' SUPSYM threshold for degeneracy and "zero" elements', 353 & ' is reset from',THRSSY,' to',THREQL 354 THRSSY = THREQL 355 END IF 356 KFRSAV = KFREE 357 ITURN = 0 358 100 ITURN = ITURN + 1 359 KFREE1 = KFREE 360C 361C Step 0: Obtain TKMO matrix for finding super symmetries 362C TKMO must be the matrix of a totally symmetric 363C operator which is very unlikely to have 364C accidental zeroes (we use the kinetic energy). 365C 366 CALL AVETK(TKMO,CMO,WRK,KFREE,LFREE) 367C Note: this can give problems for FINITE FIELD (NFIELD.gt.0) 368C because the TKMO will reflect the symmetry without the finite field(s)! 369C 370C Step 1: Find super symmetry of each MO 371C Define ISSMO(i) = super sym. of m.o. no. i 372C ISSORD(i) = sup. sym. reordering array 373C NSSYM = no. of super symmetries 374C NORBSS(*) = no. of mo's in each sup.sym. 375C IORBSS(*) = pointer array to ISSORD 376C 377 CALL IZERO(ISSMO,NORBT) 378 CALL IZERO(ISSORD,NORBT) 379 NSSYMX = 0 380 DO 180 ISYM = 1,NSYM 381 IF (NORB(ISYM) .EQ. 0) GO TO 180 382 IOFFMO = IORB(ISYM) 383 DO 170 IMO = IOFFMO+1,IOFFMO+NORB(ISYM) 384 IF (ISSMO(IMO) .EQ. 0) THEN 385 NSSYMX = NSSYMX + 1 386 IF (NSSYMX .GT. MXSSYM) GO TO 901 387C v-------- error exit 388 ISSMO(IMO) = NSSYMX 389 NORBSS(NSSYMX) = 1 390 END IF 391 ISSYM = ISSMO(IMO) 392 DO 160 JMO = IMO+1,IOFFMO+NORB(ISYM) 393 IF (ABS(TKMO(JMO,IMO)) .GT. THRSSY) THEN 394 IF (ISSMO(JMO) .EQ. 0) THEN 395C -- We have a new orbital in the supersymmetry rep. 396 ISSMO(JMO) = ISSYM 397 NORBSS(ISSYM) = NORBSS(ISSYM) + 1 398 ELSE IF (ISSMO(JMO) .NE. ISSYM) THEN 399C -- We have a conflict and must collapse the two supersym.s 400 JSSYMH = MAX(ISSMO(JMO),ISSYM) 401 JSSYML = MIN(ISSMO(JMO),ISSYM) 402 IF (IPRAVE .GE. 4) THEN 403 WRITE(LUPRI,'(/A,2I4)') 404 * ' Collapsing super symmetries :',JSSYML,JSSYMH 405 END IF 406 DO 140 KMO = IOFFMO+1,IOFFMO+NORB(ISYM) 407 IF (ISSMO(KMO) .EQ. JSSYMH) ISSMO(KMO) = JSSYML 408 140 CONTINUE 409 ISSYM = JSSYML 410 NORBSS(JSSYML) = NORBSS(JSSYML) + NORBSS(JSSYMH) 411 NORBSS(JSSYMH) = 0 412 END IF 413 END IF 414 160 CONTINUE 415 170 CONTINUE 416 180 CONTINUE 417C 418C in case of debugging, get out some print 419C 420 IF ( IPRAVE .GE. 10) THEN 421 NSSYM = NSSYMX 422 CALL IZERO(IORBSS,NSSYM) 423 WRITE(LUPRI,'(/A/A/)') 424 * ' After initializing in AVESE1 (step 1)', 425 * ' =====================================' 426 CALL AVEDMP(LUPRI) 427 END IF 428C 429C Step 2: Find NSSYM, IORBSS(*) and ISSORD(*) 430C 431 NSSYM = 0 432 IMOSS = 0 433 IORBSS(1) = 0 434 DO 260 ISSYMX = 1,NSSYMX 435 IF (NORBSS(ISSYMX) .GT. 0) THEN 436 NSSYM = NSSYM + 1 437 IF (NSSYM .NE. ISSYMX) THEN 438 NORBSS(NSSYM) = NORBSS(ISSYMX) 439 DO 220 IMO = 1,NORBT 440 IF (ISSMO(IMO) .EQ. ISSYMX) ISSMO(IMO) = NSSYM 441 220 CONTINUE 442 END IF 443 IORBSS(NSSYM) = IMOSS 444 DO 240 IMO = 1,NORBT 445 IF (ISSMO(IMO) .EQ. NSSYM) THEN 446 IMOSS = IMOSS + 1 447 ISSORD(IMOSS) = IMO 448 END IF 449 240 CONTINUE 450 IF ( (IMOSS - IORBSS(NSSYM)) .NE. NORBSS(NSSYM) ) THEN 451 GO TO 902 452 END IF 453 END IF 454 260 CONTINUE 455C 456C in case of debugging, get out some print 457C 458 IF ( IPRAVE .GE. 3) THEN 459 WRITE(LUPRI,'(/A/A/)') 460 * ' After clearing up in AVESE1 (step 2)', 461 * ' ====================================' 462 CALL AVEDMP(LUPRI) 463 END IF 464C 465C Step 3: Find degeneracies 466C 467 CALL MEMGET('INTE',KINDX1,NORBT,WRK,KFREE,LFREE) 468 CALL MEMGET('INTE',KINDX2,NORBT,WRK,KFREE,LFREE) 469 CALL AVEDEG(TKMO,WRK(KINDX1),WRK(KINDX2),NDEGSS) 470C 471C Step 4: Disentangle degenerate super symmetries 472C 473 NDIS = 0 474 DO 480 ISSYM = 1,NSSYM 475 NDEGI = NDEGSS(ISSYM) 476 IF (NDEGI .GT. 1) THEN 477 NDIS = NDIS + 1 478 IF (ITURN .EQ. 1) THEN 479 ISYM = NINFSS(ISSYM,3) 480 NORBI = NORB(ISYM) 481 NBASI = NBAS(ISYM) 482 JCMOI = ICMO(ISYM) + 1 483 NSSMOI = NORBSS(ISSYM) 484 NSSMON = NSSMOI / NDEGI 485 IF (NDEGI*NSSMON .NE. NSSMOI) THEN 486 CALL QUIT('AVESE1 error: NDEGI*NSSMON .ne. NSSMOI') 487 END IF 488 KFREE4 = KFREE 489 CALL MEMGET('REAL',KUMAT,NDEGI*NDEGI,WRK,KFREE,LFREE) 490 CALL MEMGET('INTE',KIMODE,NSSMOI,WRK,KFREE,LFREE) 491 CALL MEMGET('REAL',KWRK,NBASI*NDEGI,WRK,KFREE,LFREE) 492 CALL AVEDIS(ISSYM,ISYM,NDEGI,NORBI,NBASI,NSSMOI,NSSMON, 493 * CMO(JCMOI),TKMO,WRK(KINDX1), 494 * WRK(KUMAT),WRK(KIMODE),WRK(KWRK)) 495 CALL MEMREL('AVESE1.AVEDIS',WRK,1,KFREE4,KFREE,LFREE) 496 END IF 497 END IF 498 480 CONTINUE 499C 500 CALL MEMREL('AVESE1.step 4',WRK,1,KFREE1,KFREE,LFREE) 501C 502 IF (NDIS .GT. 0) THEN 503 IF (ITURN .EQ. 1) THEN 504 IF (IPRAVE .GE. 3) WRITE(LUPRI,'(/A)') 505 & ' AVESE1: new super symmetry analysis after'// 506 & ' disentangling degenerate super symmetries.' 507 GO TO 100 508 END IF 509 WRITE(LUPRI,'(/A)') 510 & ' AVESE1 error: disentangling failed, sorry!' 511 CALL QUIT('AVESE1 error: disentangling failed, sorry!') 512 END IF 513C 514C Step 5: Check (and change) relative phases of MO's 515C 516 IF (IPRAVE .GE. 4) THEN 517 WRITE(LUPRI,'(/A)') 518 & ' AVESE1 : checking relative phases of degenerate MOs' 519 END IF 520 IPHCHA = 0 521 DO 580 ISSYM = 1, NSSYM 522 IDEG = NINFSS(ISSYM,1) 523 ISSYMR = NINFSS(ISSYM,2) 524 IF (IDEG .GT. 1 .AND. ISSYM .EQ. ISSYMR) THEN 525C --- we have a degeneracy to other supersym.(s) 526 CALL AVEPHA(ISSYM,CMO,TKMO,IPHCHA) 527 END IF 528 580 CONTINUE 529 IF (IPRAVE .GE. 4 .AND. IPHCHA .GT. 0) THEN 530 WRITE(LUPRI,'(A,I5)') 531 & ' AVESE1 : # of phase shifts of degenerate MOs :',IPHCHA 532 END IF 533C 534C *** end of subroutine AVESE1 535C 536 CALL QEXIT('AVESE1') 537 RETURN 538C 539 901 CONTINUE 540 CALL QUIT('AVESE1: NSSYM .gt. MXSSYM in AVESET code 901') 541 902 CONTINUE 542 CALL QUIT('AVESE1: Inconsistent NORBSS(NSSYM) code 902') 543 END 544C /* Deck aveord */ 545 SUBROUTINE AVEORD() 546C 547C Aug. 2004 hjaaj - create ISSORD() from ISSMO() 548C 549C The order of the super symmetry orbitals may be changed 550C in ORDRSS and ORD2SS calls, thus we need to remake ISSORD() 551C to be sure it is correct. 552C 553#include "implicit.h" 554#include "priunit.h" 555#include "maxash.h" 556#include "maxorb.h" 557#include "infind.h" 558#include "inforb.h" 559#include "infpri.h" 560C 561 CALL QENTER('AVEORD') 562 IMOSS = 0 563 IORBSS(1) = 0 564 DO ISSYM = 1,NSSYM 565 IF (IORBSS(ISSYM) .NE. IMOSS) THEN 566 CALL AVEDMP(LUPRI) 567 CALL QUIT('AVEORD: Inconsistent IORBSS(ISSYM)') 568 END IF 569 DO IMO = 1,NORBT 570 IF (ISSMO(IMO) .EQ. ISSYM) THEN 571 IMOSS = IMOSS + 1 572 ISSORD(IMOSS) = IMO 573 END IF 574 END DO 575 IF ( (IMOSS - IORBSS(ISSYM)) .NE. NORBSS(ISSYM) ) THEN 576 CALL AVEDMP(LUPRI) 577 CALL QUIT('AVEORD: Inconsistent NORBSS(ISSYM)') 578 END IF 579 END DO 580C 581C In case of debugging, get out some print 582C 583 IF ( IPRAVE .GE. 3) THEN 584 WRITE(LUPRI,'(/A/A/)') 585 * ' After remake of ISSORD() in AVEORD ', 586 * ' ====================================' 587 CALL AVEDMP(LUPRI) 588 END IF 589 CALL QEXIT('AVEORD') 590 RETURN 591 END 592C /* Deck avetk */ 593 SUBROUTINE AVETK(TKMO,CMO,WRK,KFREE,LFREE) 594C 595C 920629-hjaaj: 596C Purpose: obtain TKMO(NORBT,NORBT) 597C 598#include "implicit.h" 599 DIMENSION CMO(*), TKMO(NORBT,NORBT), WRK(*) 600C 601 CHARACTER*8 TKLABEL 602 PARAMETER (TKLABEL = 'KINETINT') 603C 604C Used from common blocks: 605C INFORB : NORBT,NNBAST,... 606C INFPRI : IPRAVE 607C 608#include "priunit.h" 609#include "inforb.h" 610#include "infpri.h" 611C 612 LOGICAL FOUND 613C 614 CALL QENTER('AVETK ') 615 KFRSAV = KFREE 616 CALL MEMGET('REAL',KTKAO,NNBAST,WRK,KFREE,LFREE) 617 FOUND = .FALSE. 618 CALL RDONEL(TKLABEL,FOUND,WRK(KTKAO),NNBAST) 619 IF (.NOT. FOUND) THEN 620 CALL QUIT('AVETK error: TK integrals not found') 621 END IF 622C 623C 624 CALL MEMGET('REAL',KTKPK,NNORBT,WRK,KFREE,LFREE) 625 CALL MEMGET('REAL',KSCR1,NNORBX,WRK,KFREE,LFREE) 626 DO 200 ISYM = 1,NSYM 627 NORBI = NORB(ISYM) 628 IF (NORBI.EQ.0) GO TO 200 629 NBASI = NBAS(ISYM) 630 JTKAO = KTKAO + IIBAS(ISYM) 631 JTKPK = KTKPK + IIORB(ISYM) 632 JCMO = 1 + ICMO(ISYM) 633C 634 CALL UTHU(WRK(JTKAO),WRK(JTKPK),CMO(JCMO),WRK(KSCR1), 635 * NBASI,NORBI) 636C CALL UTHU(H,HT,U,WRK,NAO,NMO) 637 200 CONTINUE 638C 639 CALL PKSYM1(WRK(KSCR1),WRK(KTKPK),NORB,NSYM,-1) 640C CALL PKSYM1(ASP,APK,NDIM,NBLK,IWAY) 641C 642 CALL DSPTSI(NORBT,WRK(KSCR1),TKMO) 643C CALL DSPTSI(N,ASP,ASI) 644C 645 CALL MEMREL('AVETK',WRK,1,KFRSAV,KFREE,LFREE) 646C 647 IF ( IPRAVE .GE. 15 ) THEN 648 WRITE(LUPRI,'(3(/A))') 649 & ' Matrix used for determining supersymmetries', 650 & ' (the kinetic energy matrix)', 651 & ' ===========================================' 652 CALL OUTPUT(TKMO,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 653 END IF 654C 655 CALL QEXIT('AVETK ') 656 RETURN 657 END 658C /* Deck avedeg */ 659 SUBROUTINE AVEDEG(TKMO,INDX1,INDX2,NDEGSS) 660C 661C 920630 hjaaj 662C Purpose : find degeneracies, and construct 663C INDX1(*) = labeling degenerate orbitals 664C INDX2(*) = degeneracies of labels in INDX1 665C 666#include "implicit.h" 667 DIMENSION TKMO(NORBT,NORBT),INDX1(NORBT),INDX2(NORBT) 668 DIMENSION NDEGSS(*) 669C DIMENSION NDEGSS(MXSSYM) 670C 671C Used from common blocks: 672C INFINP : THRSSY 673C INFORB : NORBT, NSSYM, NINFSS(MXSSYM,3),MXDGSS 674C INFIND : ISSMO() 675C 676#include "maxash.h" 677#include "maxorb.h" 678#include "priunit.h" 679#include "infinp.h" 680#include "inforb.h" 681#include "infind.h" 682#include "infpri.h" 683C 684 CALL QENTER('AVEDEG') 685 CALL IZERO(INDX1,NORBT) 686 CALL IZERO(INDX2,NORBT) 687 IMOX = 0 688 NDEGMX = 0 689 DO 300 IMO = 1,NORBT 690 IF (INDX1(IMO) .EQ. 0) THEN 691 IMOX = IMOX + 1 692 INDX1(IMO) = IMOX 693 INDX2(IMOX) = 1 694 TKMOI = TKMO(IMO,IMO) 695 DO 200 JMO = IMO+1,NORBT 696 IF (ABS(TKMO(JMO,JMO)-TKMOI) .LT. THRSSY) THEN 697 INDX1(JMO) = IMOX 698 INDX2(IMOX) = INDX2(IMOX) + 1 699 END IF 700 200 CONTINUE 701 NDEGMX = MAX(NDEGMX,INDX2(IMOX)) 702 END IF 703 300 CONTINUE 704 MXDGSS = NDEGMX 705 NORXT = IMOX 706C 707C Define NINFSS(i,1) = degeneracy of MO's in supsym_i 708C NINFSS(i,2) = "root" supsym of this degeneracy 709C NINFSS(i,3) = Abelian symmetry of supsym_i 710C NDEGSS(i ) = no. of deg. MO's in this supsym 711C 712 CALL IZERO(NINFSS,MXSSYM*3) 713 CALL IZERO(NDEGSS,MXSSYM) 714 IDGERR = 0 715 DO 500 IMO = 1,NORBT 716 ISSYM = ISSMO(IMO) 717 IMOX = INDX1(IMO) 718 NDEGI = INDX2(IMOX) 719 IF (NINFSS(ISSYM,1) .EQ. 0) THEN 720 NINFSS(ISSYM,1) = NDEGI 721 NINFSS(ISSYM,2) = ISSYM 722 NINFSS(ISSYM,3) = ISMO(IMO) 723 NDEGSS(ISSYM) = 1 724 DO 400 JMO = IMO+1,NORBT 725 IF (INDX1(JMO) .EQ. IMOX) THEN 726 JSSYM = ISSMO(JMO) 727 NINFSS(JSSYM,1) = NDEGI 728 NINFSS(JSSYM,2) = ISSYM 729 NINFSS(JSSYM,3) = ISMO(JMO) 730 NDEGSS(JSSYM) = NDEGSS(JSSYM) + 1 731 END IF 732 400 CONTINUE 733 ELSE IF (NINFSS(ISSYM,1) .NE. NDEGI) THEN 734 WRITE(LUPRI,'(A,4I4)') 735 & ' Error for IMO,ISSYM, NINFSS(ISSYM,1),NDEGI = ', 736 & IMO,ISSYM, NINFSS(ISSYM,1),NDEGI 737 IDGERR = IDGERR + 1 738 END IF 739 500 CONTINUE 740C 741C Print NINFSS for debugging 742C 743 IF ( IPRAVE .GE. 3 .OR. IDGERR .GT. 0) THEN 744 IF (IDGERR .GT. 0) THEN 745 WRITE(LUPRI,'(7(/A))') 746 & ' AVEDEG: Degeneracy error, dump of NINFSS follows below.', 747 & ' Probable cause:'// 748 & ' nuclear geometry has only approximate super symmetry.', 749 & ' Some possible actions:', 750 &' (1) include more digits in nuclear coordinates', 751 &' (2) attempt to force super symmetry by making .THRSSY larger', 752 &' (3) attempt to avoid degeneracy error by making .THRSSY smaller' 753 &,' (4) disable supersymmetry by removing .SUPSYM' 754 END IF 755 WRITE(LUPRI,'(4(/A))') 756 & ' NINFSS(i,1) = degeneracy of MOs in supsym_i', 757 & ' NINFSS(i,2) = root supsym of this degeneracy', 758 & ' NINFSS(i,3) = Abelian symmetry of supsym_i', 759 & ' NDEGSS(i) = no. of deg. supsyms currently in this supsym' 760 WRITE(LUPRI,'(/A/A)') 761 & ' ISS NINFSS 1 to 3 NDEGSS', 762 & ' ============================================' 763 DO 710 ISSYM = 1, NSSYM 764 WRITE(LUPRI,'(1X,I4,4I10)') 765 & ISSYM,(NINFSS(ISSYM,J),J=1,3),NDEGSS(ISSYM) 766 710 CONTINUE 767 WRITE(LUPRI,'(A/)') 768 & ' ============================================' 769 END IF 770C 771C print for debuggers 772C 773 IF ( IPRAVE .GE. 10 .OR. IDGERR .GT. 0 ) THEN 774 WRITE(LUPRI,'(A,F15.12/)') 775 & ' Degeneracy pattern found with THRSSY =',THRSSY 776 WRITE(LUPRI,'(" I ISSMO INDX1 INDX2 TKMO(I,I)")') 777 WRITE(LUPRI,'(" ============================================")') 778 DO 700 IMO = 1, NORBT 779 WRITE(LUPRI,'(I5,3I6,F25.15)') 780 & IMO, ISSMO(IMO), INDX1(IMO), INDX2(IMO), TKMO(IMO,IMO) 781 700 CONTINUE 782 WRITE(LUPRI,'(" ============================================")') 783 END IF 784 IF (IDGERR .GT. 0) THEN 785 CALL QUIT('AVEDEG: Degeneracy error ...') 786 END IF 787C 788 CALL QEXIT('AVEDEG') 789 RETURN 790 END 791C /* Deck avedis */ 792 SUBROUTINE AVEDIS(ISSYM,ISYM,NDEGI,NORBI,NBASI,NSSMOI,NSSMON, 793 * CMOI,TKMO,INDX1,UMAT,IMODEG,WRK) 794C 795#include "implicit.h" 796C 797 DIMENSION CMOI(NBASI,NORBI), TKMO(NORBT,NORBT) 798 DIMENSION INDX1(NORBT),UMAT(NDEGI,NDEGI) 799 DIMENSION IMODEG(NDEGI,NSSMON), WRK(NBASI,NDEGI) 800C 801C Used from common blocks 802C INFORB: NSSYM, NORBT, NORBSS(), IORBSS() 803C INFIND: ISMO(),ISSORD() 804C INFPRI: IPRAVE 805C 806#include "maxash.h" 807#include "maxorb.h" 808#include "priunit.h" 809#include "inforb.h" 810#include "infind.h" 811#include "infpri.h" 812C 813 CALL QENTER('AVEDIS') 814 IF (IPRAVE .GE. 3) THEN 815 WRITE(LUPRI,'(/A,I5/A,I5)') 816 * ' AVEDIS disentangling mixed degenerate components in'// 817 * ' supersymmetry:',ISSYM,' Abelian symmetry:',ISYM 818 END IF 819 IMOXP = 0 820 KMOSS = 0 821 IMOSSE = IORBSS(ISSYM)+NORBSS(ISSYM) 822 DO 180 IMOSS = IORBSS(ISSYM)+1,IMOSSE 823 IMO = ISSORD(IMOSS) 824 IMOX = INDX1(IMO) 825 IF (IMOX .GT. IMOXP) THEN 826 IMOXP = IMOX 827 KMOSS = KMOSS + 1 828 IDEG = 1 829 IMODEG(1,KMOSS) = IMO 830 DO 120 JMOSS = IMOSS+1,IMOSSE 831 JMO = ISSORD(JMOSS) 832 IF (INDX1(JMO) .EQ. IMOX) THEN 833 IDEG = IDEG + 1 834 IMODEG(IDEG,KMOSS) = JMO 835 END IF 836 120 CONTINUE 837 END IF 838 180 CONTINUE 839 IF (KMOSS .NE. NSSMON) THEN 840 CALL QUIT('AVEDIS error: KMOSS .ne. NSSMON') 841 END IF 842C 843 IOFFMO = IORB(ISYM) 844C check that first block is diagonal 845 CALL AVEGTU(NDEGI,NORBT,1,1, 846 * IMODEG(1,1),IMODEG(1,1),UMAT,TKMO,WRK) 847 DO 280 KMOSS = 2,NSSMON 848C a) check that KMOSS block is diagonal 849 CALL AVEGTU(NDEGI,NORBT,KMOSS,KMOSS, 850 * IMODEG(1,KMOSS),IMODEG(1,KMOSS),UMAT,TKMO,WRK) 851C b) find UMAT to transform the KMOSS degenerate orbital set to 852C the same irrep rows as the first deg. orbital set 853 CALL AVEGTU(NDEGI,NORBT,1,KMOSS, 854 * IMODEG(1,1),IMODEG(1,KMOSS),UMAT,TKMO,WRK) 855 CALL DZERO(WRK,NDEGI*NBASI) 856C c) do the transformation 857 DO 240 IDEG = 1,NDEGI 858 IMO = IMODEG(IDEG,KMOSS) - IOFFMO 859 DO 220 JDEG = 1,NDEGI 860 CALL DAXPY(NBASI,UMAT(JDEG,IDEG), 861 * CMOI(1,IMO),1,WRK(1,JDEG),1) 862 220 CONTINUE 863 240 CONTINUE 864 DO 260 IDEG = 1,NDEGI 865 IMO = IMODEG(IDEG,KMOSS) - IOFFMO 866 CALL DCOPY(NBASI,WRK(1,IDEG),1,CMOI(1,IMO),1) 867 260 CONTINUE 868 280 CONTINUE 869C 870 CALL QEXIT('AVEDIS') 871 RETURN 872 END 873C /* Deck avegtu */ 874 SUBROUTINE AVEGTU(NDEGI,NORBT,IMOSS,KMOSS,IMOI,IMOK, 875 * UMAT,TKMO,UTUMAT) 876C 877#include "implicit.h" 878 DIMENSION IMOI(NDEGI), IMOK(NDEGI) 879 DIMENSION UMAT(NDEGI,NDEGI), TKMO(NORBT,NORBT) 880 DIMENSION UTUMAT(NDEGI,NDEGI) 881 PARAMETER (D1 = 1.0D0) 882C 883C Used from common blocks: 884C INFINP : THRSSY 885C INFPRI : IPRAVE 886C 887#include "maxorb.h" 888#include "priunit.h" 889#include "infinp.h" 890#include "infpri.h" 891C 892 CALL QENTER('AVEGTU') 893 DO 120 K = 1,NDEGI 894 DO 110 I = 1,NDEGI 895 UMAT(I,K) = TKMO( IMOI(I) , IMOK(K) ) 896 110 CONTINUE 897 120 CONTINUE 898C 899C -- Do some printout for debuggers 900C 901 IF (IPRAVE .GE. 10) THEN 902 WRITE(LUPRI,'(/A/A/A,3I5)') 903 * ' Test output from AVEGTU', 904 * ' =======================', 905 * ' IMOSS, KMOSS, NDEGI =',IMOSS,KMOSS,NDEGI 906 WRITE(LUPRI,'(A,20I4)') ' IMOI :',(IMOI(I),I=1,NDEGI) 907 WRITE(LUPRI,'(A,20I4)') ' IMOK :',(IMOK(I),I=1,NDEGI) 908 WRITE(LUPRI,'(/A/A/)') ' UMAT found in AVEGTU', 909 * ' ====================' 910 CALL OUTPUT(UMAT,1,NDEGI,1,NDEGI,NDEGI,NDEGI,1,LUPRI) 911 END IF 912C 913 IF (IMOSS .EQ. KMOSS) THEN 914C ... check that UMAT is diagonal 915 NERR = NOFDIA(NDEGI,NDEGI,UMAT,THRSSY) 916 IF (NERR .GT. 0) THEN 917 WRITE(LUPRI,'(2A,I5)') ' AVEGTU error, ', 918 * 'no. of non-zero off-diagonal elements in U =',NERR 919 WRITE(LUPRI,'(A,I4)') ' -- IMOSS = ', IMOSS 920 WRITE(LUPRI,'(A,I4)') ' -- KMOSS = ', KMOSS 921 WRITE(LUPRI,'(/A)') ' The U matrix :' 922 CALL OUTPUT(UMAT,1,NDEGI,1,NDEGI,NDEGI,NDEGI,1,LUPRI) 923 CALL QUIT('AVEGTU error for IMOSS .eq. KMOSS') 924 END IF 925 ELSE 926C ... IMOSS .ne. KMOSS 927C a) check that matrix is x*orthogonal matrix 928 ITURN = 0 929 300 ITURN = ITURN + 1 930 CALL DGEMM('T','N',NDEGI,NDEGI,NDEGI,1.D0, 931 & UMAT,NDEGI, 932 & UMAT,NDEGI,0.D0, 933 & UTUMAT,NDEGI) 934 NERR = NOFDIA(NDEGI,NDEGI,UTUMAT,THRSSY) 935 IF (NERR .GT. 0 .OR. IPRAVE .GE. 10) THEN 936 IF (NERR .GT. 0) WRITE(LUPRI,'(2A,I5)') ' AVEGTU error, ', 937 * 'no. of non-zero off-diagonal elements in U^(T)U =',NERR 938 WRITE(LUPRI,'(A,I4)') ' -- IMOSS = ', IMOSS 939 WRITE(LUPRI,'(A,I4)') ' -- KMOSS = ', KMOSS 940 WRITE(LUPRI,'(/A)') ' The U matrix :' 941 CALL OUTPUT(UMAT,1,NDEGI,1,NDEGI,NDEGI,NDEGI,1,LUPRI) 942 WRITE(LUPRI,'(/A)') ' The UTU matrix :' 943 CALL OUTPUT(UTUMAT,1,NDEGI,1,NDEGI,NDEGI,NDEGI,1,LUPRI) 944 IF (NERR .GT. 0) 945 & CALL QUIT('AVEGTU error for IMOSS .ne. KMOSS') 946 END IF 947C 948C HJMAERKE 940818-hjaaj: we ought to force UMAT orthogonal with 949C symmetric orthonormalization; with the current code we may 950C get deviations from orthonormality up to THRSSY 951 X2 = UTUMAT(1,1) 952 DO 310 I = 2,NDEGI 953 X2 = X2 * UTUMAT(I,I) 954 310 CONTINUE 955C det(U) = sqrt(X2) = X2 ** ( .5 ) 956 XEXP = 2*NDEGI 957 XEXP = - D1 / XEXP 958 XINV = X2 ** ( XEXP ) 959 CALL DSCAL(NDEGI*NDEGI,XINV,UMAT,1) 960 IF (ABS(D1-XINV) .GT. THRSSY) THEN 961 IF (IPRAVE .GE. 10 .OR. ITURN .GE. 3) THEN 962 WRITE(LUPRI,'(A,I4)') ' -- ITURN = ', ITURN 963 WRITE(LUPRI,'(A,1D20.10)') ' -- XINV = ', XINV 964 WRITE(LUPRI,'(A,I4)') ' -- IMOSS = ', IMOSS 965 WRITE(LUPRI,'(A,I4)') ' -- KMOSS = ', KMOSS 966 WRITE(LUPRI,'(/A)') ' The U matrix :' 967 CALL OUTPUT(UMAT,1,NDEGI,1,NDEGI,NDEGI,NDEGI,1,LUPRI) 968 END IF 969 IF (ITURN .LT. 3) GO TO 300 970 WRITE(LUPRI,'(A)') 971 * ' AVEGTU error : failed to normalize UMAT' 972 CALL QUIT('AVEGTU error: could not normalize UMAT') 973 END IF 974 END IF 975 CALL QEXIT('AVEGTU') 976 RETURN 977 END 978C /* Deck avepha */ 979 SUBROUTINE AVEPHA(ISSYMR,CMO,TKMO,IPHCHA) 980C 981C 2-Jul-1992 hjaaj+hh 982C Purpose: make relative phases the same in each block 983C of a set of degenerate mo.s. 984C A block corresponds to orbitals following 985C a row of an irrep. 986C 987#include "implicit.h" 988 DIMENSION CMO(NCMOT), TKMO(NORBT,NORBT) 989C 990C Used from common blocks: 991C INFINP : THRSSY 992C INFORB : NCMOT,NORBT,NORBSS(), IORBSS(),...,NINFSS(*,3),... 993C INFIND : ISSORD() 994C 995#include "maxash.h" 996#include "maxorb.h" 997#include "priunit.h" 998#include "infinp.h" 999#include "inforb.h" 1000#include "infind.h" 1001#include "infpri.h" 1002C 1003 CALL QENTER('AVEPHA') 1004 IF (NINFSS(ISSYMR,2) .NE. ISSYMR) THEN 1005 WRITE(LUPRI,'(/A)') 1006 * ' AVEPHA error: ISSYMR is not a "root" symmetry' 1007 WRITE(LUPRI,'(A,I4)') ' ISSYMR = ', ISSYMR 1008 WRITE(LUPRI,'(A,I4)') ' NINFSS(ISSYMR,2) = ', NINFSS(ISSYMR,2) 1009 CALL QUIT('AVEPHA error: ISSYMR is not a "root" symmetry') 1010 END IF 1011 IOSSI = IORBSS(ISSYMR) 1012 NOSSI = NORBSS(ISSYMR) 1013 DO 800 JSSYM = 1,NSSYM 1014 JSSYMR = NINFSS(JSSYM,2) 1015 IF (JSSYMR .EQ. ISSYMR) THEN 1016 IF (JSSYM .LT. ISSYMR) THEN 1017 CALL QUIT('AVEPHA error: JSSYM .lt ISSYMR') 1018 END IF 1019 IF (NORBSS(JSSYM) .NE. NOSSI) THEN 1020 WRITE(LUPRI,'(/A/A,2I5/A)') 1021 * ' AVEHPA error: Diff. no. of orb.s in deg. supsym.s', 1022 * ' JSSYMR, JSSYM =',JSSYMR,JSSYM, 1023 * ' Dump of super symmetry information :' 1024 CALL AVEDMP(LUPRI) 1025 CALL QUIT( 1026 & 'AVEHPA error: Diff. no. of orb.s in deg. supsym.s') 1027 END IF 1028 ISYMJ = NINFSS(JSSYM,3) 1029 ICMOJ = ICMO(ISYMJ) 1030 NBASJ = NBAS(ISYMJ) 1031 IORBJ = IORB(ISYMJ) 1032 IOSSJ = IORBSS(JSSYM) 1033 IMO1 = ISSORD(IOSSI + 1) 1034 JMO1 = ISSORD(IOSSJ + 1) 1035 DO 700 K = 1,NOSSI 1036 IMO = ISSORD(IOSSI + K) 1037 JMO = ISSORD(IOSSJ + K) 1038 IF (ABS(TKMO(IMO,IMO) - TKMO(JMO,JMO)) .GT. THRSSY) THEN 1039 WRITE (LUPRI,'(/A,2I5)') 1040 & ' AVEPHA error: different ordering'// 1041 & ' of degenerate orbitals in supsyms:',ISSYMR,JSSYM 1042 WRITE(LUPRI,'(/A,2I5/A,D20.10/A,D20.10)') 1043 & ' IMO, JMO =',IMO,JMO, 1044 & ' TKMO(IMO,IMO) =',TKMO(IMO,IMO), 1045 & ' TKMO(JMO,JMO) =',TKMO(JMO,JMO) 1046 CALL QUIT( 1047 & 'AVEPHA error: different ordering of deg. orb.s') 1048 END IF 1049 IF (ABS(TKMO(IMO,IMO1)) .LT. THRSSY) THEN 1050 WRITE(LUPRI,'(/A/A,2I5/A,D10.2)') 1051 & ' AVEPHA error: TKMO(IMO,IMO1) is zero', 1052 & ' IMO, IMO1 =',IMO,IMO1, 1053 & ' TKMO(IMO,IMO1) =',TKMO(IMO,IMO1) 1054 CALL QTRACE(LUERR) 1055 CALL QUIT('AVEPHA error: TKMO(IMO,IMO1) is zero') 1056 ELSE IF (ABS(TKMO(IMO,IMO1)+TKMO(JMO,JMO1)) .LT. THRSSY) 1057 * THEN 1058 IPHCHA = IPHCHA + 1 1059 IF (IPRAVE .GE. 4) THEN 1060 WRITE(LUPRI,'(A,I4)') 1061 * ' AVEPHA changing phase of mo no.',JMO 1062 END IF 1063 JOFF = ICMOJ + NBASJ*(JMO-IORBJ-1) 1064 DO 600 J = 1,NBASJ 1065 CMO(JOFF+J) = -CMO(JOFF+J) 1066 600 CONTINUE 1067 ELSE IF (ABS(TKMO(IMO,IMO1)-TKMO(JMO,JMO1)) .GT. THRSSY) 1068 * THEN 1069 WRITE(LUPRI,'(/A,2(/A,2I5/A,D20.10))') 1070 & ' AVEPHA error: TKMO(IMO,IMO1) not TKMO(JMO,JMO1)', 1071 & ' IMO, IMO1 =',IMO,IMO1, 1072 & ' TKMO(IMO,IMO1) =',TKMO(IMO,IMO1), 1073 & ' JMO, JMO1 =',JMO,JMO1, 1074 & ' TKMO(JMO,JMO1) =',TKMO(JMO,JMO1) 1075 CALL QTRACE(LUERR) 1076 CALL QUIT('AVEPHA error: error in deg. orb.s code 2') 1077 END IF 1078 700 CONTINUE 1079 END IF 1080 800 CONTINUE 1081 CALL QEXIT('AVEPHA') 1082 RETURN 1083 END 1084#if defined (VAR_NOTUSED ) 1085 AVEUDV is not used 1086C /* Deck aveudv */ 1087 SUBROUTINE AVEUDV(UDV) 1088C 1089C 09-Jul-1992 Hinne Hettema 1090C 1091C Purpose : perform an averaging of the one-electron density 1092C matrix. First we identify an IMO which is degenerate with 1093C other MO's, then we average the matrix rows and colums. 1094C 1095#include "implicit.h" 1096 DIMENSION UDV(NASHT,*) 1097C 1098C Used from common blocks 1099C INFORB : NSYM ,NSSYM, NORBSS(*), IORBSS(*), NINFSS(*,3) 1100C INFIND : ICH(*), ISSMO(*), ISSORD(*) 1101C 1102#include "maxorb.h" 1103#include "maxash.h" 1104#include "inforb.h" 1105#include "infind.h" 1106#include "infpri.h" 1107C 1108 PARAMETER ( MAXDEG = 20, D0 = 0.0D0, D1 = 1.0D0 ) 1109 DIMENSION ISSOFF(MAXDEG),IDEGV(MAXDEG,2) 1110C 1111 CALL QENTER('AVEUDV') 1112C 1113 DO 800 ISSYMR = 1, NSSYM 1114C -- skip if supersym not degenerate or not 'root' supersym 1115 IF ( NINFSS(ISSYMR,1) .EQ. 1) GO TO 800 1116 IF ( NINFSS(ISSYMR,2) .NE. ISSYMR) GO TO 800 1117C 1118 NDEG = 0 1119 DO 100 JSSYM = ISSYMR, NSSYM 1120C -- skip if not same 'root' supersymmetry 1121C (because then not degenerate at all) 1122 IF ( NINFSS(JSSYM,2) .EQ. ISSYMR) THEN 1123 NDEG = NDEG + 1 1124 ISSOFF(NDEG) = IORBSS(JSSYM) - IORBSS(ISSYMR) 1125 END IF 1126 100 CONTINUE 1127 IF ( NDEG .NE. NINFSS(ISSYMR,1) ) THEN 1128 WRITE(LUERR,'(A)') ' AVEUDV error : inconsistent'// 1129 & ' degeneracies' 1130 WRITE(LUERR,'(A,I4)') ' - ISSYMR = ', ISSYMR 1131 WRITE(LUERR,'(A,I4)') ' - NDEG = ', NDEG 1132 WRITE(LUERR,'(A,I4)') ' - NINFSS(ISSYMR,1) = ', 1133 & NINFSS(ISSYMR,1) 1134 CALL QUIT(' AVEUDV : inconsistent degeneracy in NINFSS') 1135 END IF 1136C 1137C -- Compute averaging factor 1138 AVFAC = NDEG 1139 AVFAC = D1 / AVFAC 1140C 1141C -- in the following, use supersymmetry info to get IMO 1142 IRSYM = NINFSS(ISSYMR,3) 1143 NASHI = NASH(IRSYM) 1144 IOFFSR = IORBSS(ISSYMR) 1145 IOFFDV = IASH(IRSYM) + 1 1146 DO 700 IMOSS = IOFFSR+1, IOFFSR+NORBSS(ISSYMR) 1147 IMO = ISSORD(IMOSS) 1148 IMOW = ICH(IMO) 1149 IF (IMOW .LE. 0) GO TO 700 1150 DO 600 JMOSS = IMOSS, IOFFSR+NORBSS(ISSYMR) 1151 JMO = ISSORD(JMOSS) 1152 JMOW = ICH(JMO) 1153 IF ( JMOW .LE. 0 ) GO TO 600 1154 DO 300 IDEG = 2, NDEG 1155 KMO = ISSORD(IMOSS + ISSOFF(IDEG)) 1156 KMOW = ICH(KMO) 1157 IF (KMOW .LE. 0 ) THEN 1158 CALL QUIT('AVEUDV code 300.1: ICH(KMO).le.0') 1159 END IF 1160 LMO = ISSORD(JMOSS + ISSOFF(IDEG)) 1161 LMOW = ICH(LMO) 1162 IF (LMOW .LE. 0 ) THEN 1163 CALL QUIT('AVEUDV code 300.2: ICH(LMO).le.0') 1164 END IF 1165 IDEGV(IDEG,1) = KMOW 1166 IDEGV(IDEG,2) = LMOW 1167 300 CONTINUE 1168C -- now collect all elements which should be made equal 1169 XUDVIJ = UDV(IMOW,JMOW) 1170 DO 400 IDEG = 2, NDEG 1171 KMOW = IDEGV(IDEG,1) 1172 LMOW = IDEGV(IDEG,2) 1173 XUDVIJ = XUDVIJ + UDV(KMOW,LMOW) 1174 400 CONTINUE 1175 XUDVIJ = XUDVIJ * AVFAC 1176C -- now distribute all elements which should be made equal 1177 UDV(IMOW,JMOW) = XUDVIJ 1178 IF (IMOW .NE. JMOW) UDV(JMOW,IMOW) = XUDVIJ 1179 DO 500 IDEG = 2, NDEG 1180 KMOW = IDEGV(IDEG,1) 1181 LMOW = IDEGV(IDEG,2) 1182 UDV(KMOW,LMOW) = XUDVIJ 1183 IF(KMOW .NE. LMOW) UDV(LMOW,KMOW) = XUDVIJ 1184 500 CONTINUE 1185 600 CONTINUE 1186 700 CONTINUE 1187C 1188 800 CONTINUE 1189C 1190C Zero elements which are not totally symmetric 1191C with respect to the supersymmetry. 1192C 1193 DO 980 ISYM = 1,NSYM 1194 IASHI = IASH(ISYM) 1195 NASHI = NASH(ISYM) 1196 DO 960 IMOW = IASHI+1,IASHI+NASHI 1197 IMO = ISX(NISHT+IMOW) 1198 DO 950 JMOW = IMOW+1,IASHI+NASHI 1199 JMO = ISX(NISHT+JMOW) 1200 IF (ISSMO(JMO) .NE. ISSMO(IMO)) THEN 1201 UDV(IMOW,JMOW) = D0 1202 UDV(JMOW,IMOW) = D0 1203 END IF 1204 950 CONTINUE 1205 960 CONTINUE 1206 980 CONTINUE 1207C 1208 CALL QEXIT('AVEUDV') 1209C 1210 RETURN 1211 END 1212C /* Deck aveh1 */ 1213 AVEH1 is not used 1214 SUBROUTINE AVEH1(H1) 1215C 1216C 12-Jul-1992 Hinne Hettema 1217C 1218C Purpose : perform an averaging of a one-electron 1219C matrix. First we identify an IMO which is degenerate with 1220C other MO's, then we average the matrix rows and colums. 1221C 1222#include "implicit.h" 1223 DIMENSION H1(NORBT,*) 1224C 1225C Used from common blocks 1226C INFORB : NSYM ,NSSYM, NORBSS(*), IORBSS(*), NINFSS(*,3) 1227C INFIND : ICH(*), ISSMO(*), ISSORD(*) 1228C 1229#include "maxorb.h" 1230#include "maxash.h" 1231#include "priunit.h" 1232#include "inforb.h" 1233#include "infind.h" 1234#include "infpri.h" 1235C 1236 PARAMETER ( MAXDEG = 20, D0 = 0.0D0, D1 = 1.0D0 ) 1237 DIMENSION ISSOFF(MAXDEG),IDEGV(MAXDEG,2) 1238C 1239 CALL QENTER('AVEH1 ') 1240C 1241 IF ( IPRAVE .GE. 10 ) THEN 1242 WRITE(LUPRI,'(//A)') ' H1 matrix before averaging' 1243 WRITE(LUPRI,'(A)') ' ==========================' 1244 CALL OUTPUT(H1,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 1245 END IF 1246C 1247 DO 800 ISSYMR = 1, NSSYM 1248C -- skip if supersym not degenerate or not 'root' supersym 1249 IF ( NINFSS(ISSYMR,1) .EQ. 1) GO TO 800 1250 IF ( NINFSS(ISSYMR,2) .NE. ISSYMR) GO TO 800 1251C 1252 NDEG = 0 1253 DO 100 JSSYM = ISSYMR, NSSYM 1254C -- skip if not same 'root' supersymmetry 1255C (because then not degenerate at all) 1256 IF ( NINFSS(JSSYM,2) .EQ. ISSYMR) THEN 1257 NDEG = NDEG + 1 1258 ISSOFF(NDEG) = IORBSS(JSSYM) - IORBSS(ISSYMR) 1259 END IF 1260 100 CONTINUE 1261 IF ( NDEG .NE. NINFSS(ISSYMR,1) ) THEN 1262 WRITE(LUERR,'(A)') ' AVEH1 error : inconsistent'// 1263 & ' degeneracies' 1264 WRITE(LUERR,'(A,I4)') ' - ISSYMR = ', ISSYMR 1265 WRITE(LUERR,'(A,I4)') ' - NDEG = ', NDEG 1266 WRITE(LUERR,'(A,I4)') ' - NINFSS(ISSYMR,1) = ', 1267 & NINFSS(ISSYMR,1) 1268 CALL QUIT(' AVEH1 : inconsistent degeneracy in NINFSS') 1269 END IF 1270C 1271C -- Compute averaging factor 1272 AVFAC = NDEG 1273 AVFAC = D1 / AVFAC 1274C 1275C -- in the following, use supersymmetry info to get IMO 1276 IRSYM = NINFSS(ISSYMR,3) 1277 NORBI = NORB(IRSYM) 1278 IOFFSR = IORBSS(ISSYMR) 1279 DO 700 IMOSS = IOFFSR+1, IOFFSR+NORBSS(ISSYMR) 1280 IMO = ISSORD(IMOSS) 1281 DO 600 JMOSS = IMOSS, IOFFSR+NORBSS(ISSYMR) 1282 JMO = ISSORD(JMOSS) 1283 IDEGV(1,1) = IMO 1284 IDEGV(1,2) = JMO 1285 DO 300 IDEG = 2, NDEG 1286 IDEGV(IDEG,1) = ISSORD(IMOSS + ISSOFF(IDEG)) 1287 IDEGV(IDEG,2) = ISSORD(JMOSS + ISSOFF(IDEG)) 1288 300 CONTINUE 1289C -- now collect all elements which should be made equal 1290 XH1IJ = H1(IMO,JMO) 1291 DO 400 IDEG = 2, NDEG 1292 KMO = IDEGV(IDEG,1) 1293 LMO = IDEGV(IDEG,2) 1294 XH1IJ = XH1IJ + H1(KMO,LMO) 1295 400 CONTINUE 1296 XH1IJ = XH1IJ * AVFAC 1297C -- now distribute all elements which should be made equal 1298 H1(IMO,JMO) = XH1IJ 1299 IF (IMO .NE. JMO) H1(JMO,IMO) = XH1IJ 1300 DO 500 IDEG = 2, NDEG 1301 KMO = IDEGV(IDEG,1) 1302 LMO = IDEGV(IDEG,2) 1303 H1(KMO,LMO) = XH1IJ 1304 IF(KMO .NE. LMO) H1(LMO,KMO) = XH1IJ 1305 500 CONTINUE 1306 600 CONTINUE 1307 700 CONTINUE 1308C 1309 800 CONTINUE 1310C 1311C Zero elements which are not totally symmetric 1312C with respect to the supersymmetry. 1313C 1314 DO 980 ISYM = 1,NSYM 1315 IORBI = IORB(ISYM) 1316 NORBI = NORB(ISYM) 1317 DO 960 IMO = IORBI+1,IORBI+NORBI 1318 DO 950 JMO = IMO+1,IORBI+NORBI 1319 IF (ISSMO(JMO) .NE. ISSMO(IMO)) THEN 1320 H1(IMO,JMO) = D0 1321 H1(JMO,IMO) = D0 1322 END IF 1323 950 CONTINUE 1324 960 CONTINUE 1325 980 CONTINUE 1326C 1327 IF ( IPRAVE .GE. 10 ) THEN 1328 WRITE(LUPRI,'(//A)') ' H1 matrix after averaging' 1329 WRITE(LUPRI,'(A)') ' =========================' 1330 CALL OUTPUT(H1,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 1331 END IF 1332C 1333 CALL QEXIT('AVEH1 ') 1334C 1335 RETURN 1336 END 1337#endif 1338C /* Deck avedv */ 1339 SUBROUTINE AVEDV(DV) 1340C 1341C 09-Jul-1992 Hinne Hettema 1342C 1343C Purpose : perform an averaging of the one-electron density 1344C matrix. First we identify an IMO which is degenerate with 1345C other MO's, then we average the matrix rows and colums. 1346C 1347#include "implicit.h" 1348 DIMENSION DV(NNASHX) 1349C 1350C Used from common blocks 1351C INFORB : NSYM ,NSSYM, NORBSS(*), IORBSS(*), NINFSS(*,3) 1352C INFIND : IROW(*), ICH(*), ISSMO(*), ISSORD(*) 1353C 1354#include "maxorb.h" 1355#include "maxash.h" 1356#include "priunit.h" 1357#include "inforb.h" 1358#include "infind.h" 1359#include "infpri.h" 1360C 1361 PARAMETER ( MAXDEG = 20, D0 = 0.0D0, D1 = 1.0D0 ) 1362 DIMENSION ISSOFF(MAXDEG),IDEGV(MAXDEG) 1363C 1364 CALL QENTER('AVEDV ') 1365C 1366 IF (IPRAVE .GE. 5) THEN 1367 WRITE(LUPRI,'(/A)') 1368 * ' AVEDV test output: DV matrix before averaging' 1369 CALL OUTPAK(DV,NASHT,1,LUPRI) 1370 END IF 1371C 1372 DO 800 ISSYMR = 1, NSSYM 1373C -- skip if supersym not degenerate or not 'root' supersym 1374 IF ( NINFSS(ISSYMR,1) .EQ. 1) GO TO 800 1375 IF ( NINFSS(ISSYMR,2) .NE. ISSYMR) GO TO 800 1376C 1377 NDEG = 0 1378 DO 100 JSSYM = ISSYMR, NSSYM 1379C -- skip if not same 'root' supersymmetry 1380C (because then not degenerate at all) 1381 IF ( NINFSS(JSSYM,2) .EQ. ISSYMR) THEN 1382 NDEG = NDEG + 1 1383 ISSOFF(NDEG) = IORBSS(JSSYM) - IORBSS(ISSYMR) 1384 END IF 1385 100 CONTINUE 1386 IF ( NDEG .NE. NINFSS(ISSYMR,1) ) THEN 1387 WRITE(LUERR,'(A)') ' AVEDV error : inconsistent'// 1388 & ' degeneracies' 1389 WRITE(LUERR,'(A,I4)') ' - ISSYMR = ', ISSYMR 1390 WRITE(LUERR,'(A,I4)') ' - NDEG = ', NDEG 1391 WRITE(LUERR,'(A,I4)') ' - NINFSS(ISSYMR,1) = ', 1392 & NINFSS(ISSYMR,1) 1393 CALL QUIT(' AVEDV : inconsistent degeneracy in NINFSS') 1394 END IF 1395C 1396C -- Compute averaging factor 1397 AVFAC = NDEG 1398 AVFAC = D1 / AVFAC 1399C 1400C -- in the following, use supersymmetry info to get IMO 1401 IRSYM = NINFSS(ISSYMR,3) 1402 NASHI = NASH(IRSYM) 1403 IOFFSR = IORBSS(ISSYMR) 1404 IOFFDV = IASH(IRSYM) + 1 1405 DO 700 IMOSS = IOFFSR+1, IOFFSR+NORBSS(ISSYMR) 1406 IMO = ISSORD(IMOSS) 1407 IMOW = ICH(IMO) 1408 IF (IMOW .LE. 0) GO TO 700 1409 DO 600 JMOSS = IOFFSR+1, IMOSS 1410 JMO = ISSORD(JMOSS) 1411 JMOW = ICH(JMO) 1412 IF ( JMOW .LE. 0 ) GO TO 600 1413 IDEGV(1) = IROW(IMOW) + JMOW 1414 DO 300 IDEG = 2, NDEG 1415 KMO = ISSORD(IMOSS + ISSOFF(IDEG)) 1416 KMOW = ICH(KMO) 1417 IF (KMOW .LE. 0 ) THEN 1418 CALL QUIT('AVEDV code 300.1') 1419 END IF 1420 LMO = ISSORD(JMOSS + ISSOFF(IDEG)) 1421 LMOW = ICH(LMO) 1422 IF (LMOW .LE. 0 ) THEN 1423 CALL QUIT('AVEDV code 300.2') 1424 END IF 1425 IDEGV(IDEG) = IROW(KMOW) + LMOW 1426 300 CONTINUE 1427C -- now collect all elements which should be made equal 1428 XDVIJ = DV(IDEGV(1)) 1429 DO 400 IDEG = 2, NDEG 1430 XDVIJ = XDVIJ + DV(IDEGV(IDEG)) 1431 400 CONTINUE 1432 XDVIJ = XDVIJ * AVFAC 1433C -- now distribute all elements which should be made equal 1434 DO 500 IDEG = 1, NDEG 1435 DV(IDEGV(IDEG)) = XDVIJ 1436 500 CONTINUE 1437 600 CONTINUE 1438 700 CONTINUE 1439C 1440 800 CONTINUE 1441C 1442C Zero elements which are not totally symmetric 1443C with respect to the supersymmetry. 1444C 1445 DO 980 ISYM = 1,NSYM 1446 IASHI = IASH(ISYM) 1447 NASHI = NASH(ISYM) 1448 DO 960 IMOW = IASHI+2,IASHI+NASHI 1449 IMO = ISX(NISHT+IMOW) 1450 DO 950 JMOW = IASHI+1,IMOW-1 1451 JMO = ISX(NISHT+JMOW) 1452 IF (ISSMO(JMO) .NE. ISSMO(IMO)) THEN 1453 IDV = IROW(IMOW) + JMOW 1454 DV(IDV) = D0 1455 END IF 1456 950 CONTINUE 1457 960 CONTINUE 1458 980 CONTINUE 1459C 1460 IF (IPRAVE .GE. 5) THEN 1461 WRITE(LUPRI,'(/A)') 1462 * ' AVEDV test output: DV matrix after averaging' 1463 CALL OUTPAK(DV,NASHT,1,LUPRI) 1464 END IF 1465 CALL QEXIT('AVEDV ') 1466C 1467 RETURN 1468 END 1469#if defined (VAR_NOTUSED ) 1470C /* Deck aveh1p */ 1471 SUBROUTINE AVEH1P(H1P) 1472C 1473C 12-Jul-1992 Hinne Hettema 1474C 1475C Purpose : perform an averaging of a one-electron 1476C matrix. First we identify an IMO which is degenerate with 1477C other MO's, then we average the matrix rows and colums. 1478C 1479#include "implicit.h" 1480 DIMENSION H1P(NNORBX) 1481C 1482C Used from common blocks 1483C INFORB : NSYM ,NSSYM, NORBSS(*), IORBSS(*), NINFSS(*,3) 1484C INFIND : ISSMO(*), ISSORD(*) 1485C 1486#include "maxorb.h" 1487#include "maxash.h" 1488#include "priunit.h" 1489#include "inforb.h" 1490#include "infind.h" 1491#include "infpri.h" 1492C 1493 PARAMETER ( MAXDEG = 20, D0 = 0.0D0, D1 = 1.0D0 ) 1494 DIMENSION ISSOFF(MAXDEG),IDEGV(MAXDEG) 1495C 1496 CALL QENTER('AVEH1P') 1497C 1498 IF (IPRAVE .GE. 5) THEN 1499 WRITE(LUPRI,'(/A)') 1500 * ' AVEH1P test output: H1P matrix before averaging' 1501 CALL OUTPAK(H1P,NORBT,LUPRI) 1502 END IF 1503C 1504 DO 800 ISSYMR = 1, NSSYM 1505C -- skip if supersym not degenerate or not 'root' supersym 1506 IF ( NINFSS(ISSYMR,1) .EQ. 1) GO TO 800 1507 IF ( NINFSS(ISSYMR,2) .NE. ISSYMR) GO TO 800 1508C 1509 NDEG = 0 1510 DO 100 JSSYM = ISSYMR, NSSYM 1511C -- skip if not same 'root' supersymmetry 1512C (because then not degenerate at all) 1513 IF ( NINFSS(JSSYM,2) .EQ. ISSYMR) THEN 1514 NDEG = NDEG + 1 1515 ISSOFF(NDEG) = IORBSS(JSSYM) - IORBSS(ISSYMR) 1516 END IF 1517 100 CONTINUE 1518 IF ( NDEG .NE. NINFSS(ISSYMR,1) ) THEN 1519 WRITE(LUERR,'(A)') ' AVEH1P error : inconsistent'// 1520 & ' degeneracies' 1521 WRITE(LUERR,'(A,I4)') ' - ISSYMR = ', ISSYMR 1522 WRITE(LUERR,'(A,I4)') ' - NDEG = ', NDEG 1523 WRITE(LUERR,'(A,I4)') ' - NINFSS(ISSYMR,1) = ', 1524 & NINFSS(ISSYMR,1) 1525 CALL QUIT(' AVEH1P : inconsistent degeneracy in NINFSS') 1526 END IF 1527C 1528C -- Compute averaging factor 1529 AVFAC = NDEG 1530 AVFAC = D1 / AVFAC 1531C 1532C -- in the following, use supersymmetry info to get IMO 1533 IRSYM = NINFSS(ISSYMR,3) 1534 NORBI = NORB(IRSYM) 1535 IOFFSR = IORBSS(ISSYMR) 1536 DO 700 IMOSS = IOFFSR+1, IOFFSR+NORBSS(ISSYMR) 1537 IMO = ISSORD(IMOSS) 1538 DO 600 JMOSS = IOFFSR+1, IMOSS 1539 JMO = ISSORD(JMOSS) 1540 IDEGV(1) = IROW(IMO) + JMO 1541 DO 300 IDEG = 2, NDEG 1542C 920714-hjaaj: KMO var. is used to be safe 1543C xlf version 2.2 has made errors in a similar statem. 1544C to the version with IROW(IMO+ISSOFF(IDEG)) 1545 KMO = ISSORD(IMOSS + ISSOFF(IDEG)) 1546 IDEGV(IDEG) = IROW(KMO) + ISSORD(JMOSS + ISSOFF(IDEG)) 1547 300 CONTINUE 1548C -- now collect all elements which should be made equal 1549 XH1PIJ = H1P(IDEGV(1)) 1550 DO 400 IDEG = 2, NDEG 1551 XH1PIJ = XH1PIJ + H1P(IDEGV(IDEG)) 1552 400 CONTINUE 1553 XH1PIJ = XH1PIJ * AVFAC 1554C -- now distribute all elements which should be made equal 1555 DO 500 IDEG = 1, NDEG 1556 H1P(IDEGV(IDEG)) = XH1PIJ 1557 500 CONTINUE 1558 600 CONTINUE 1559 700 CONTINUE 1560C 1561 800 CONTINUE 1562C 1563C Zero elements which are not totally symmetric 1564C with respect to the supersymmetry. 1565C 1566 DO 980 ISYM = 1,NSYM 1567 IORBI = IORB(ISYM) 1568 NORBI = NORB(ISYM) 1569 DO 960 IMO = IORBI+2,IORBI+NORBI 1570 DO 950 JMO = IORBI+1,IMO-1 1571 IF (ISSMO(JMO) .NE. ISSMO(IMO)) THEN 1572 IJH1P = IROW(IMO) + JMO 1573 H1P(IJH1P) = D0 1574 END IF 1575 950 CONTINUE 1576 960 CONTINUE 1577 980 CONTINUE 1578C 1579 IF (IPRAVE .GE. 5) THEN 1580 WRITE(LUPRI,'(/A)') 1581 * ' AVEH1P test output: H1P matrix after averaging' 1582 CALL OUTPAK(H1P,NORBT,1,LUPRI) 1583 END IF 1584C 1585 CALL QEXIT('AVEH1P') 1586C 1587 RETURN 1588 END 1589#endif 1590C /* Deck avefck */ 1591 SUBROUTINE AVEFCK(FCK) 1592C 1593C 12-Jul-1992 Hinne Hettema 1594C 1595C Purpose : perform an averaging of a one-electron 1596C matrix. First we identify an IMO which is degenerate with 1597C other MO's, then we average the matrix rows and colums. 1598C 1599#include "implicit.h" 1600 DIMENSION FCK(NNORBT) 1601C 1602C Used from common blocks 1603C INFORB : NSYM ,NSSYM, NORBSS(*), IORBSS(*), NINFSS(*,3) 1604C INFIND : ISSMO(*), ISSORD(*) 1605C 1606#include "maxorb.h" 1607#include "maxash.h" 1608#include "priunit.h" 1609#include "inforb.h" 1610#include "infind.h" 1611#include "infpri.h" 1612C 1613 PARAMETER ( MAXDEG = 20, D0 = 0.0D0, D1 = 1.0D0 ) 1614 DIMENSION ISMDEG(MAXDEG),ISSOFF(MAXDEG),IDEGV(MAXDEG) 1615C 1616 CALL QENTER('AVEFCK') 1617C 1618 IF (IPRAVE .GE. 7) THEN 1619 WRITE(LUPRI,'(/A)') 1620 * ' AVEFCK test output: FCK matrix before averaging' 1621 CALL OUTPKB(FCK,NORB,NSYM,1,LUPRI) 1622 END IF 1623C 1624 DO 800 ISSYMR = 1, NSSYM 1625C -- skip if supersym not degenerate or not 'root' supersym 1626 IF ( NINFSS(ISSYMR,1) .EQ. 1) GO TO 800 1627 IF ( NINFSS(ISSYMR,2) .NE. ISSYMR) GO TO 800 1628C 1629 NDEG = 0 1630 DO 100 JSSYM = ISSYMR, NSSYM 1631C -- skip if not same 'root' supersymmetry 1632C (because then not degenerate at all) 1633 IF ( NINFSS(JSSYM,2) .EQ. ISSYMR) THEN 1634 NDEG = NDEG + 1 1635 ISMDEG(NDEG) = NINFSS(JSSYM,3) 1636 ISSOFF(NDEG) = IORBSS(JSSYM) - IORBSS(ISSYMR) 1637 END IF 1638 100 CONTINUE 1639 IF ( NDEG .NE. NINFSS(ISSYMR,1) ) THEN 1640 WRITE(LUERR,'(A)') ' AVEFCK error : inconsistent'// 1641 & ' degeneracies' 1642 WRITE(LUERR,'(A,I4)') ' - ISSYMR = ', ISSYMR 1643 WRITE(LUERR,'(A,I4)') ' - NDEG = ', NDEG 1644 WRITE(LUERR,'(A,I4)') ' - NINFSS(ISSYMR,1) = ', 1645 & NINFSS(ISSYMR,1) 1646 CALL QUIT(' AVEFCK : inconsistent degeneracy in NINFSS') 1647 END IF 1648C 1649C -- Compute averaging factor 1650 AVFAC = NDEG 1651 AVFAC = D1 / AVFAC 1652C 1653C -- in the following, use supersymmetry info to get IMO 1654 IRSYM = ISMDEG(1) 1655 IORBI = IORB(IRSYM) 1656 NORBI = NORB(IRSYM) 1657 IIORBI = IIORB(IRSYM) 1658 IOFFSR = IORBSS(ISSYMR) 1659 DO 700 IMOSS = IOFFSR+1, IOFFSR+NORBSS(ISSYMR) 1660 IMO = ISSORD(IMOSS) - IORBI 1661 DO 600 JMOSS = IOFFSR+1, IMOSS 1662 JMO = ISSORD(JMOSS) - IORBI 1663 IDEGV(1) = IIORBI + IROW(IMO) + JMO 1664 DO 300 IDEG = 2, NDEG 1665 JRSYM = ISMDEG(IDEG) 1666 IORBJ = IORB(JRSYM) 1667 IIORBJ = IIORB(JRSYM) 1668C 920714-hjaaj: KMO var. is used to be safe 1669C xlf version 2.2 has made errors in a similar statem. 1670C to the version with IROW(IMO+ISSOFF(IDEG)) 1671 KMO = ISSORD(IMOSS + ISSOFF(IDEG)) - IORBJ 1672 LMO = ISSORD(JMOSS + ISSOFF(IDEG)) - IORBJ 1673 IDEGV(IDEG) = IIORBJ + IROW(KMO) + LMO 1674 300 CONTINUE 1675C -- now collect all elements which should be made equal 1676 XFCKIJ = FCK(IDEGV(1)) 1677 DO 400 IDEG = 2, NDEG 1678 XFCKIJ = XFCKIJ + FCK(IDEGV(IDEG)) 1679 400 CONTINUE 1680 XFCKIJ = XFCKIJ * AVFAC 1681C -- now distribute all elements which should be made equal 1682 DO 500 IDEG = 1, NDEG 1683 FCK(IDEGV(IDEG)) = XFCKIJ 1684 500 CONTINUE 1685 600 CONTINUE 1686 700 CONTINUE 1687C 1688 800 CONTINUE 1689C 1690C Zero elements which are not totally symmetric 1691C with respect to the supersymmetry. 1692C 1693 DO 980 ISYM = 1,NSYM 1694 IORBI = IORB(ISYM) 1695 NORBI = NORB(ISYM) 1696 IIORBI = IIORB(ISYM) 1697 DO 960 IMO = 2,NORBI 1698 DO 950 JMO = 1,IMO-1 1699 IF (ISSMO(IORBI+JMO) .NE. ISSMO(IORBI+IMO)) THEN 1700 IJFCK = IIORBI + IROW(IMO) + JMO 1701 FCK(IJFCK) = D0 1702 END IF 1703 950 CONTINUE 1704 960 CONTINUE 1705 980 CONTINUE 1706C 1707 IF (IPRAVE .GE. 7) THEN 1708 WRITE(LUPRI,'(/A)') 1709 * ' AVEFCK test output: FCK matrix after averaging' 1710 CALL OUTPKB(FCK,NORB,NSYM,1,LUPRI) 1711 END IF 1712C 1713 CALL QEXIT('AVEFCK') 1714C 1715 RETURN 1716 END 1717C /* Deck aveprt */ 1718 SUBROUTINE AVEPRT(LUPRI) 1719C 1720C 7-Jul-92 Hinne Hettema and Hans Joergen Aa. Jensen. 1721C Revised 15-Jul-1993 hjaaj 1722C Purpose : write supersymmetry information 1723C 1724#include "implicit.h" 1725C 1726#include "maxash.h" 1727#include "maxorb.h" 1728#include "inforb.h" 1729#include "infind.h" 1730C 1731 PARAMETER ( ISTEP = 10 ) 1732C 1733 CALL QENTER('AVEPRT') 1734C 1735C Check if supersymmetry is identical to "D2h" symmetry 1736C 1737 IF (NSSYM .LE. NSYM .AND. MXDGSS .EQ. 1) THEN 1738 DO 200 ISYM = 1,NSYM 1739 JSYM = 0 1740 DO 100 ISSYM = 1,NSSYM 1741 IF (NINFSS(ISSYM,3) .EQ. ISYM) JSYM = JSYM + 1 1742 100 CONTINUE 1743 IF (JSYM .GT. 1) GO TO 300 1744 200 CONTINUE 1745 WRITE(LUPRI,'(/5X,A)') 'No supersymmetry found.' 1746 GO TO 9999 1747 END IF 1748 300 CONTINUE 1749C 1750 WRITE(LUPRI,'(/5X,A)') 'Supersymmetry specification' 1751 WRITE(LUPRI,'( 5X,A)') '===========================' 1752 IF (MXDGSS .EQ. 1) 1753 &WRITE(LUPRI,'( 5X,A)') '(only one dimensional irreps)' 1754 ICOUNT = 0 1755 DO 500 ISSYM = 1, NSSYM, ISTEP 1756 ICOUNT = ICOUNT + 1 1757 ISTRT = ISSYM 1758 IEND = MIN(NSSYM, ICOUNT * ISTEP) 1759 ILEN = IEND - ISTRT + 1 1760C 1761 WRITE(LUPRI,1100) 1762 & 'Supersymmetry number', 1763 & (I, I=ISTRT, IEND) 1764 IF (MXDGSS .GT. 1) WRITE(LUPRI,1100) 1765 & 'degenerate with sup.sym.', 1766 & (NINFSS(I,2), I = ISTRT, IEND) 1767 WRITE(LUPRI,1100) 1768 & 'and located in abelian symmetry ', 1769 & (NINFSS(I,3), I=ISTRT, IEND) 1770 WRITE(LUPRI,1110) 1771 & (' --', I = 1, ILEN) 1772 WRITE(LUPRI,1100) 1773 & 'Total number of orbitals', 1774 & (NORBSS(I), I = ISTRT, IEND) 1775 WRITE(LUPRI,'()') 1776 500 CONTINUE 1777 9999 CALL QEXIT('AVEPRT') 1778 RETURN 1779 1100 FORMAT(5X,A,T36,10I4) 1780 1110 FORMAT(5X, T36,10A4) 1781 END 1782C /* Deck avecph */ 1783 SUBROUTINE AVECPH(IPHCHA,CMO,WRK,KFRSAV,LFRSAV) 1784C 1785C 920715-hjaaj+hh 1786C 1787C Purpose: change phases such that degenerate MO's have 1788C the same relative phase. 1789C 1790#include "implicit.h" 1791 DIMENSION CMO(NCMOT), WRK(*) 1792C 1793C Used from common blocks: 1794C INFORB : NCMOT,N2ORBX,NINFSS(*,4), NSSYM 1795C 1796#include "inforb.h" 1797C 1798 CALL QENTER('AVECPH') 1799 KFREE = KFRSAV 1800 LFREE = LFRSAV 1801 CALL MEMGET('REAL',KTKMO,N2ORBX,WRK,KFREE,LFREE) 1802 CALL AVETK(WRK(KTKMO),CMO,WRK,KFREE,LFREE) 1803 IPHCHA = 0 1804 DO 100 ISSYM = 1, NSSYM 1805 IDEG = NINFSS(ISSYM,1) 1806 ISSYMR = NINFSS(ISSYM,2) 1807 IF (IDEG .GT. 1 .AND. ISSYM .EQ. ISSYMR) THEN 1808C --- we have a degeneracy to other supersym.(s) 1809 CALL AVEPHA(ISSYMR,CMO,WRK(KTKMO),IPHCHA) 1810 END IF 1811 100 CONTINUE 1812 CALL MEMREL('AVECPH',WRK,1,KFRSAV,KFREE,LFREE) 1813 CALL QEXIT('AVECPH') 1814 RETURN 1815 END 1816C /* Deck avechk */ 1817 SUBROUTINE AVECHK(INPERR) 1818C 1819C 14-Jul-1992 Hinne Hettema and Hans Joergen Aa. Jensen. 1820C 1821C Purpose: Check supersymmetry information which has been built up 1822C up till now. 1823C 1824#include "implicit.h" 1825C 1826C Used from common blocks: 1827C INFINP : SUPSYM 1828C INFORB : NSSYM,NORBSS(),IORBSS(),NINFSS(*,3) 1829C INFIND : ISW(*),ISSORD(*) 1830C 1831#include "maxash.h" 1832#include "maxorb.h" 1833#include "priunit.h" 1834#include "infinp.h" 1835#include "infind.h" 1836#include "inforb.h" 1837#include "infpri.h" 1838#include "orbtypdef.h" 1839C 1840 PARAMETER (MAXDEG = 20, D1 = 1.0D0) 1841 DIMENSION ISSDEG(MAXDEG),ISSOFF(MAXDEG) 1842C 1843 CALL QENTER('AVECHK') 1844 IF (.NOT.SUPSYM) GO TO 9999 1845C 1846 DO 800 ISSYMR = 1,NSSYM 1847C Skip this supsym if not degenerate or not "root supsym" 1848 IF (NINFSS(ISSYMR,1) .EQ. 1) GO TO 800 1849 IF (NINFSS(ISSYMR,2) .NE. ISSYMR) GO TO 800 1850 NDEG = 0 1851 DO 100 ISSYM = ISSYMR,NSSYM 1852 IF (NINFSS(ISSYM,2) .EQ. ISSYMR) THEN 1853 NDEG = NDEG + 1 1854 ISSDEG(NDEG) = ISSYM 1855 ISSOFF(NDEG) = IORBSS(ISSYM) - IORBSS(ISSYMR) 1856 END IF 1857 100 CONTINUE 1858 IF (NDEG .NE. NINFSS(ISSYMR,1)) THEN 1859 INPERR = INPERR + 1 1860 WRITE(LUPRI,'(/A/A,I4/A,I4,A,I4)') 1861 & ' Super symmetry error: inconsistent degeneracy in NINFSS', 1862 & ' for "root" super symmetry',ISSYMR, 1863 & ' Found:',NDEG,'; Expected:',NINFSS(ISSYMR,1) 1864 END IF 1865C 1866 IOFFSR = IORBSS(ISSYMR) 1867 ITYPA = 0 1868 DO 480 IMOSS = IOFFSR+1,IOFFSR+NORBSS(ISSYMR) 1869 IMO = ISSORD(IMOSS) 1870 ITYP = IOBTYP(IMO) 1871 IF (ITYP .EQ. JTACT) THEN 1872 ITYPA = IACTYP( ICH(IMO) ) 1873 END IF 1874 DO 410 IDEG = 2,NDEG 1875 KMO = ISSORD(IMOSS + ISSOFF(IDEG)) 1876 KTYP = IOBTYP(KMO) 1877 IF ( ITYP .NE. KTYP ) THEN 1878 INPERR = INPERR + 1 1879 WRITE (LUPRI,'(/A,2I5/A,2X,A,2X,A/A,I4)') 1880 & ' Super symmetry error: degenerate orbitals',IMO,KMO, 1881 & ' are not of same type',COBTYP(ITYP),COBTYP(KTYP), 1882 & ' for "root" super symmetry',ISSYMR 1883 ELSE IF (ITYP .EQ. JTACT) THEN 1884 KTYPA = IACTYP( ICH(KMO) ) 1885 IF (ITYPA .NE. KTYPA) THEN 1886 INPERR = INPERR + 1 1887 WRITE (LUPRI,'(/A,2I5/A,2I5/A,I4)') 1888 & ' Super symmetry error: degenerate active orbitals',IMO,KMO, 1889 & ' are not in same RAS shell type',ITYPA,KTYPA, 1890 & ' for "root" super symmetry',ISSYMR 1891 END IF 1892 END IF 1893 410 CONTINUE 1894 480 CONTINUE 1895 800 CONTINUE 1896C 1897 9999 CONTINUE 1898 CALL QEXIT('AVECHK') 1899 RETURN 1900 END 1901C /* Deck avedmp */ 1902 SUBROUTINE AVEDMP(LUWDMP) 1903C 1904C 2-Jul-92 Hinne Hettema 1905C Purpose : dump supersymmetry information to output for debugging 1906C 1907#include "implicit.h" 1908C 1909#include "maxash.h" 1910#include "maxorb.h" 1911#include "inforb.h" 1912#include "infind.h" 1913C 1914 CALL QENTER('AVEDMP') 1915 WRITE(LUWDMP,'(/A/A)') 1916 * ' Dump of super symmetry information', 1917 * ' ========================================' 1918 WRITE(LUWDMP,'(/A,I4)') 1919 * ' Number of supersymmetries (NSSYM) is ',NSSYM 1920 WRITE(LUWDMP,'(/A/A)') 1921 * ' ISS NORBSS IORBSS', 1922 * ' ==============================' 1923 DO 100 I = 1, NSSYM 1924 WRITE(LUWDMP,'(1X,I4,2I10)') I, NORBSS(I), IORBSS(I) 1925 100 CONTINUE 1926 WRITE(LUWDMP,'(A)') ' ==============================' 1927 WRITE(LUWDMP,'(/A/A)') 1928 * ' I ISSMO ISSORD ISMO', 1929 * ' ========================================' 1930 DO 200 I = 1, NORBT 1931 WRITE(LUWDMP,'(1X,I4,3I10)') I, ISSMO(I), ISSORD(I), ISMO(I) 1932 200 CONTINUE 1933 WRITE(LUWDMP,'(A/)') ' ========================================' 1934 CALL QEXIT('AVEDMP') 1935 RETURN 1936 END 1937