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 19c /* deck CC_3HYP */ 20*=====================================================================* 21 SUBROUTINE CC_3HYP(WORK,LWORK) 22*---------------------------------------------------------------------* 23* 24* Purpose: direct calculation of frequency-dependent third 25* hyperpolarizabilities for the Couple Cluster models 26* 27* CCS, CC2, CCSD 28* 29* Written by Christof Haettig, april 1997. 30* 31*=====================================================================* 32#if defined (IMPLICIT_NONE) 33 IMPLICIT NONE 34#else 35# include "implicit.h" 36#endif 37#include "priunit.h" 38#include "ccsdinp.h" 39#include "cc4rinf.h" 40#include "ccr1rsp.h" 41 42* local parameters: 43 CHARACTER*(19) MSGDBG 44 PARAMETER (MSGDBG = '[debug] CC_3HYP> ') 45 LOGICAL LOCDBG 46 PARAMETER (LOCDBG = .FALSE. ) 47 48 CHARACTER*10 MODEL 49 LOGICAL LADD 50 INTEGER LWORK 51 INTEGER KHYPPOL, NBHYPPOL 52 INTEGER K0HTRAN, K0HDOTS, K0HCONS, N0HTRAN, MX0HTRAN, MX0HDOTS 53 INTEGER K1HTRAN, K1HDOTS, K1HCONS, N1HTRAN, MX1HTRAN, MX1HDOTS 54 INTEGER K0GTRAN, K0GDOTS, K0GCONS, N0GTRAN, MX0GTRAN, MX0GDOTS 55 INTEGER K1GTRAN, K1GDOTS, K1GCONS, N1GTRAN, MX1GTRAN, MX1GDOTS 56 INTEGER K2GTRAN, K2GDOTS, K2GCONS, N2GTRAN, MX2GTRAN, MX2GDOTS 57 INTEGER K1FTRAN, K1FDOTS, K1FCONS, N1FTRAN, MX1FTRAN, MX1FDOTS 58 INTEGER K2FTRAN, K2FDOTS, K2FCONS, N2FTRAN, MX2FTRAN, MX2FDOTS 59 INTEGER K0FATRAN,K0FADOTS,K0FACONS,N0FATRAN,MX0FATRAN,MX0FADOTS 60 INTEGER K1FATRAN,K1FADOTS,K1FACONS,N1FATRAN,MX1FATRAN,MX1FADOTS 61 INTEGER K2FATRAN,K2FADOTS,K2FACONS,N2FATRAN,MX2FATRAN,MX2FADOTS 62 INTEGER K2EATRAN,K2EADOTS,K2EACONS,N2EATRAN,MX2EATRAN,MX2EADOTS 63 INTEGER KXTRAN, KXDOTS, KXCONS, NXTRAN, MXXTRAN, MXXDOTS 64 INTEGER MXTRAN3, MXTRAN4, MXVEC1, MXVEC2 65 INTEGER KEND0, LEND0, IOPT 66 67#if defined (SYS_CRAY) 68 REAL WORK(LWORK) 69#else 70 DOUBLE PRECISION WORK(LWORK) 71#endif 72 73*---------------------------------------------------------------------* 74* print header for hyperpolarizability section 75*---------------------------------------------------------------------* 76 WRITE (LUPRI,'(7(/1X,2A),/)') 77 & '************************************', 78 & '*******************************', 79 & '* ', 80 & ' *', 81 & '*---------- OUTPUT FROM COUPLED CLUSTE', 82 & 'R QUARTIC RESPONSE ---------*', 83 & '* ', 84 & ' *', 85 & '*-------- CALCULATION OF CC HYPER', 86 & 'POLARIZABILITIES ---------*', 87 & '* ', 88 & ' *', 89 & '************************************', 90 & '*******************************' 91 92*---------------------------------------------------------------------* 93 IF (.NOT. (CCS .OR. CC2 .OR. CCSD) ) THEN 94 CALL QUIT('CC_3HYP called for unknown Coupled Cluster.') 95 END IF 96 97* print some debug/info output 98 IF (IPR4HYP .GT. 10) WRITE(LUPRI,*) 'CC_3HYP Workspace:',LWORK 99 100*---------------------------------------------------------------------* 101* allocate & initialize work space for hyperpolarizabilities 102*---------------------------------------------------------------------* 103 NBHYPPOL = N4ROPER * N4RFREQ 104 105 KHYPPOL = 1 106 KEND0 = KHYPPOL + 2*NBHYPPOL 107 108 109 MXTRAN3 = 2 * NLRTLBL * NLRTLBL * NLRTLBL 110 MXVEC2 = (NLRTLBL+1) * NLRTLBL / 2 111 MXTRAN4 = 2 * NLRTLBL * NLRTLBL * NLRTLBL * NLRTLBL 112 MXVEC1 = NLRTLBL 113 MXTRAN3 = MIN(MXTRAN3,2*60*NBHYPPOL) ! 60 = # perm. for F^A{O} 114 MXTRAN4 = MIN(MXTRAN4,2*30*NBHYPPOL) ! 30 = # perm. for F^AB{O} 115 MXVEC2 = MIN(MXVEC2,2*60*NBHYPPOL) ! 60 = # perm. for F^A{O} 116 MXVEC1 = MIN(MXVEC1,2*30*NBHYPPOL) ! 30 = # perm. for F^AB{O} 117 118 MX0HTRAN = 5 * MXTRAN3 119 MX1HTRAN = 5 * MXTRAN4 120 MX0GTRAN = 4 * MXTRAN3 121 MX1GTRAN = 4 * MXTRAN3 122 MX2GTRAN = 4 * MXTRAN4 123 MX1FTRAN = 3 * MXTRAN3 124 MX2FTRAN = 3 * MXTRAN3 125 MX0FATRAN = 5 * MXTRAN3 126 MX1FATRAN = 5 * MXTRAN3 127 MX2FATRAN = 5 * MXTRAN4 128 MX2EATRAN = 3 * MXTRAN3 129 MXXTRAN = 1 * MXTRAN3 130 131 MX0HDOTS = MXTRAN3 * MXVEC2 132 MX1HDOTS = MXTRAN4 * MXVEC1 133 MX0GDOTS = MXTRAN3 * MXVEC2 134 MX1GDOTS = MXTRAN3 * MXVEC2 135 MX2GDOTS = MXTRAN4 * MXVEC1 136 MX1FDOTS = MXTRAN3 * MXVEC2 137 MX2FDOTS = MXTRAN3 * MXVEC2 138 MX0FADOTS = MXTRAN3 * MXVEC2 139 MX1FADOTS = MXTRAN3 * MXVEC2 140 MX2FADOTS = MXTRAN4 * MXVEC1 141 MX2EADOTS = MXTRAN3 * MXVEC2 142 MXXDOTS = MXTRAN3 * MXVEC2 143 144 K0HTRAN = KEND0 145 K0HDOTS = K0HTRAN + MX0HTRAN 146 K0HCONS = K0HDOTS + MX0HDOTS 147 KEND0 = K0HCONS + MX0HDOTS 148 149 K1HTRAN = KEND0 150 K1HDOTS = K1HTRAN + MX1HTRAN 151 K1HCONS = K1HDOTS + MX1HDOTS 152 KEND0 = K1HCONS + MX1HDOTS 153 154 K0GTRAN = KEND0 155 K0GDOTS = K0GTRAN + MX0GTRAN 156 K0GCONS = K0GDOTS + MX0GDOTS 157 KEND0 = K0GCONS + MX0GDOTS 158 159 K1GTRAN = KEND0 160 K1GDOTS = K1GTRAN + MX1GTRAN 161 K1GCONS = K1GDOTS + MX1GDOTS 162 KEND0 = K1GCONS + MX1GDOTS 163 164 K2GTRAN = KEND0 165 K2GDOTS = K2GTRAN + MX2GTRAN 166 K2GCONS = K2GDOTS + MX2GDOTS 167 KEND0 = K2GCONS + MX2GDOTS 168 169 K1FTRAN = KEND0 170 K1FDOTS = K1FTRAN + MX1FTRAN 171 K1FCONS = K1FDOTS + MX1FDOTS 172 KEND0 = K1FCONS + MX1FDOTS 173 174 K2FTRAN = KEND0 175 K2FDOTS = K2FTRAN + MX2FTRAN 176 K2FCONS = K2FDOTS + MX2FDOTS 177 KEND0 = K2FCONS + MX2FDOTS 178 179 K0FATRAN = KEND0 180 K0FADOTS = K0FATRAN + MX0FATRAN 181 K0FACONS = K0FADOTS + MX0FADOTS 182 KEND0 = K0FACONS + MX0FADOTS 183 184 K1FATRAN = KEND0 185 K1FADOTS = K1FATRAN + MX1FATRAN 186 K1FACONS = K1FADOTS + MX1FADOTS 187 KEND0 = K1FACONS + MX1FADOTS 188 189 K2FATRAN = KEND0 190 K2FADOTS = K2FATRAN + MX2FATRAN 191 K2FACONS = K2FADOTS + MX2FADOTS 192 KEND0 = K2FACONS + MX2FADOTS 193 194 K2EATRAN = KEND0 195 K2EADOTS = K2EATRAN + MX2EATRAN 196 K2EACONS = K2EADOTS + MX2EADOTS 197 KEND0 = K2EACONS + MX2EADOTS 198 199 KXTRAN = KEND0 200 KXDOTS = KXTRAN + MXXTRAN 201 KXCONS = KXDOTS + MXXDOTS 202 KEND0 = KXCONS + MXXDOTS 203 204 LEND0 = LWORK - KEND0 205 206 IF (LEND0.LT.0) THEN 207 WRITE (LUPRI,*) 'KEND0,LEND0:',KEND0,LEND0 208 CALL QUIT('Insufficient memory in CC_3HYP.') 209 END IF 210 211* initialize hyperpolarizabilities and the lists of contributions: 212 CALL DZERO(WORK(KHYPPOL),2*NBHYPPOL) 213 CALL DZERO(WORK(K0HTRAN),MX0HTRAN+2*MX0HDOTS) 214 CALL DZERO(WORK(K1HTRAN),MX1HTRAN+2*MX1HDOTS) 215 CALL DZERO(WORK(K0GTRAN),MX0GTRAN+2*MX0GDOTS) 216 CALL DZERO(WORK(K1GTRAN),MX1GTRAN+2*MX1GDOTS) 217 CALL DZERO(WORK(K2GTRAN),MX2GTRAN+2*MX2GDOTS) 218 CALL DZERO(WORK(K1FTRAN),MX1FTRAN+2*MX1FDOTS) 219 CALL DZERO(WORK(K2FTRAN),MX2FTRAN+2*MX2FDOTS) 220 CALL DZERO(WORK(K0FATRAN),MX0FATRAN+2*MX0FADOTS) 221 CALL DZERO(WORK(K1FATRAN),MX1FATRAN+2*MX1FADOTS) 222 CALL DZERO(WORK(K2FATRAN),MX2FATRAN+2*MX2FADOTS) 223 CALL DZERO(WORK(K2EATRAN),MX2EATRAN+2*MX2EADOTS) 224 CALL DZERO(WORK(KXTRAN),MXXTRAN+2*MXXDOTS) 225 226*---------------------------------------------------------------------* 227* set up lists for H, G and F transformations and ETA vectors: 228*---------------------------------------------------------------------* 229 230 LADD = .FALSE. 231 232 CALL CC4R_SETUP(MXTRAN3, MXVEC2, MXTRAN4, MXVEC1, 233 & WORK(K0HTRAN), WORK(K0HDOTS), WORK(K0HCONS), N0HTRAN, 234 & WORK(K1HTRAN), WORK(K1HDOTS), WORK(K1HCONS), N1HTRAN, 235 & WORK(K0GTRAN), WORK(K0GDOTS), WORK(K0GCONS), N0GTRAN, 236 & WORK(K1GTRAN), WORK(K1GDOTS), WORK(K1GCONS), N1GTRAN, 237 & WORK(K2GTRAN), WORK(K2GDOTS), WORK(K2GCONS), N2GTRAN, 238 & WORK(K1FTRAN), WORK(K1FDOTS), WORK(K1FCONS), N1FTRAN, 239 & WORK(K2FTRAN), WORK(K2FDOTS), WORK(K2FCONS), N2FTRAN, 240 & WORK(K0FATRAN),WORK(K0FADOTS),WORK(K0FACONS),N0FATRAN, 241 & WORK(K1FATRAN),WORK(K1FADOTS),WORK(K1FACONS),N1FATRAN, 242 & WORK(K2FATRAN),WORK(K2FADOTS),WORK(K2FACONS),N2FATRAN, 243 & WORK(K2EATRAN),WORK(K2EADOTS),WORK(K2EACONS),N2EATRAN, 244 & WORK(KXTRAN), WORK(KXDOTS), WORK(KXCONS), NXTRAN, 245 & WORK(KHYPPOL), NBHYPPOL, LADD ) 246 247*---------------------------------------------------------------------* 248* calculate H matrix contributions: 249*---------------------------------------------------------------------* 250 IF (.NOT.L_USE_CHI3) THEN 251 IOPT = 5 252 CALL CC_HMAT('L0','R1','R1','R1','R2',N0HTRAN, MXVEC2, 253 & WORK(K0HTRAN),WORK(K0HDOTS),WORK(K0HCONS), 254 & WORK(KEND0), LEND0, IOPT ) 255 END IF 256 257 IOPT = 5 258 CALL CC_HMAT('L1 ','R1 ','R1 ','R1 ','R1 ',N1HTRAN, MXVEC1, 259 & WORK(K1HTRAN),WORK(K1HDOTS),WORK(K1HCONS), 260 & WORK(KEND0), LEND0, IOPT ) 261 262*---------------------------------------------------------------------* 263* calculate G matrix contributions: 264*---------------------------------------------------------------------* 265 IOPT = 5 266 CALL CC_GMATRIX('L0 ','R2 ','R1 ','R2 ',N0GTRAN, MXVEC2, 267 & WORK(K0GTRAN),WORK(K0GDOTS),WORK(K0GCONS), 268 & WORK(KEND0), LEND0, IOPT ) 269 270 IF (L_USE_CHI3) CALL DSCAL(MX0GDOTS,-1.0d0,WORK(K0GCONS),1) 271 272 IF (.NOT.L_USE_CHI3) THEN 273 IOPT = 5 274 CALL CC_GMATRIX('L1 ','R1 ','R1 ','R2 ',N1GTRAN, MXVEC2, 275 & WORK(K1GTRAN),WORK(K1GDOTS),WORK(K1GCONS), 276 & WORK(KEND0), LEND0, IOPT ) 277 END IF 278 279 IOPT = 5 280 CALL CC_GMATRIX('L2 ','R1 ','R1 ','R1 ',N2GTRAN, MXVEC1, 281 & WORK(K2GTRAN),WORK(K2GDOTS),WORK(K2GCONS), 282 & WORK(KEND0), LEND0, IOPT ) 283 284 285*---------------------------------------------------------------------* 286* calculate F matrix contributions: 287*---------------------------------------------------------------------* 288 IOPT = 5 289 CALL CC_FMATRIX(WORK(K1FTRAN),N1FTRAN,'L1 ','R2 ',IOPT,'R2 ', 290 & WORK(K1FDOTS),WORK(K1FCONS),MXVEC2, 291 & WORK(KEND0), LEND0) 292 293 IF (L_USE_CHI3) CALL DSCAL(MX1FDOTS,-1.0d0,WORK(K1FCONS),1) 294 295 IF (.NOT.L_USE_CHI3) THEN 296 IOPT = 5 297 CALL CC_FMATRIX(WORK(K2FTRAN),N2FTRAN,'L2 ','R1 ',IOPT,'R2 ', 298 & WORK(K2FDOTS),WORK(K2FCONS),MXVEC2, 299 & WORK(KEND0), LEND0) 300 END IF 301 302*---------------------------------------------------------------------* 303* calculate F{O} matrix contributions: 304*---------------------------------------------------------------------* 305 CALL CCQR_FADRV('L0 ','o1 ','R2 ','R2 ',N0FATRAN, MXVEC2, 306 & WORK(K0FATRAN),WORK(K0FADOTS), 307 & WORK(K0FACONS), 308 & WORK(KEND0), LEND0, 'DOTP' ) 309 310 IF (L_USE_CHI3) CALL DSCAL(MX0FADOTS,-1.0d0,WORK(K0FACONS),1) 311 312 IF (.NOT.L_USE_CHI3) THEN 313 CALL CCQR_FADRV('L1 ','o1 ','R1 ','R2 ',N1FATRAN, MXVEC2, 314 & WORK(K1FATRAN),WORK(K1FADOTS), 315 & WORK(K1FACONS), 316 & WORK(KEND0), LEND0, 'DOTP' ) 317 END IF 318 319 CALL CCQR_FADRV('L2 ','o1 ','R1 ','R1 ',N2FATRAN, MXVEC1, 320 & WORK(K2FATRAN),WORK(K2FADOTS), 321 & WORK(K2FACONS), 322 & WORK(KEND0), LEND0, 'DOTP' ) 323 324*---------------------------------------------------------------------* 325* calculate ETA{O} vector contributions: 326*---------------------------------------------------------------------* 327 IF (.NOT.L_USE_CHI3) THEN 328 CALL CCQR_EADRV('L2 ','o1 ','R2 ',N2EATRAN, MXVEC2, 329 & WORK(K2EATRAN),WORK(K2EADOTS),WORK(K2EACONS), 330 & WORK(KEND0), LEND0, 'DOTP' ) 331 END IF 332 333*---------------------------------------------------------------------* 334* chi vector contributions: 335*---------------------------------------------------------------------* 336 IF (L_USE_CHI3) THEN 337 CALL CC_DOTDRV('X3 ','R2 ',NXTRAN,MXVEC2, 338 & WORK(KXTRAN), WORK(KXDOTS), WORK(KXCONS), 339 & WORK(KEND0), LEND0 ) 340C CALL CC_DOTDRV('L3 ','O2 ',NXTRAN,MXVEC2, 341C & WORK(KXTRAN), WORK(KXDOTS), WORK(KXCONS), 342C & WORK(KEND0), LEND0 ) 343 END IF 344 345*---------------------------------------------------------------------* 346* collect contributions and add them to hyperpolarizabilities 347*---------------------------------------------------------------------* 348 349 LADD = .TRUE. 350 351 CALL CC4R_SETUP(MXTRAN3, MXVEC2, MXTRAN4, MXVEC1, 352 & WORK(K0HTRAN), WORK(K0HDOTS), WORK(K0HCONS), N0HTRAN, 353 & WORK(K1HTRAN), WORK(K1HDOTS), WORK(K1HCONS), N1HTRAN, 354 & WORK(K0GTRAN), WORK(K0GDOTS), WORK(K0GCONS), N0GTRAN, 355 & WORK(K1GTRAN), WORK(K1GDOTS), WORK(K1GCONS), N1GTRAN, 356 & WORK(K2GTRAN), WORK(K2GDOTS), WORK(K2GCONS), N2GTRAN, 357 & WORK(K1FTRAN), WORK(K1FDOTS), WORK(K1FCONS), N1FTRAN, 358 & WORK(K2FTRAN), WORK(K2FDOTS), WORK(K2FCONS), N2FTRAN, 359 & WORK(K0FATRAN),WORK(K0FADOTS),WORK(K0FACONS),N0FATRAN, 360 & WORK(K1FATRAN),WORK(K1FADOTS),WORK(K1FACONS),N1FATRAN, 361 & WORK(K2FATRAN),WORK(K2FADOTS),WORK(K2FACONS),N2FATRAN, 362 & WORK(K2EATRAN),WORK(K2EADOTS),WORK(K2EACONS),N2EATRAN, 363 & WORK(KXTRAN), WORK(KXDOTS), WORK(KXCONS), NXTRAN, 364 & WORK(KHYPPOL), NBHYPPOL, LADD ) 365 366 367*---------------------------------------------------------------------* 368* print output & return: 369*---------------------------------------------------------------------* 370 371 CALL CC4HYPPRT(WORK(KHYPPOL)) 372 373 RETURN 374 END 375 376*=====================================================================* 377* END OF SUBROUTINE CC_3HYP * 378*=====================================================================* 379c /* deck CC4HYPPRT */ 380*=====================================================================* 381 SUBROUTINE CC4HYPPRT(HYPERPOL) 382*---------------------------------------------------------------------* 383* 384* Purpose: print third hyperpolarizabilities 385* 386* Written by Christof Haettig in april 1997. 387* 388*=====================================================================* 389#if defined (IMPLICIT_NONE) 390 IMPLICIT NONE 391#else 392# include "implicit.h" 393#endif 394#include "priunit.h" 395#include "ccorb.h" 396#include "ccsdinp.h" 397#include "cc4rinf.h" 398#include "ccroper.h" 399 400 401 CHARACTER*10 BLANKS 402 CHARACTER*80 STRING 403 INTEGER ISYMA, ISYMB, ISYMC, ISYMD, ISYME 404 INTEGER IFREQ, IOPER, ISYOLD, ISYTOT 405 406 407#if defined (SYS_CRAY) 408 REAL HYPERPOL(N4RFREQ,N4ROPER,2) 409 REAL HALF 410#else 411 DOUBLE PRECISION HYPERPOL(N4RFREQ,N4ROPER,2) 412 DOUBLE PRECISION HALF 413#endif 414 PARAMETER (HALF = 0.5D0) 415 416*---------------------------------------------------------------------* 417* print header for hyperpolarizability section 418*---------------------------------------------------------------------* 419 BLANKS = ' ' 420 STRING = ' RESULTS FOR THE THIRD HYPERPOLARIZABILITIES ' 421 422 IF (CCS) THEN 423 CALL AROUND( BLANKS//'FINAL CCS'//STRING(1:45)//BLANKS ) 424 ELSE IF (CC2) THEN 425 CALL AROUND( BLANKS//'FINAL CC2'//STRING(1:45)//BLANKS ) 426 ELSE IF (CCSD) THEN 427 CALL AROUND( BLANKS//'FINAL CCSD'//STRING(1:45)//BLANKS ) 428 ELSE 429 CALL QUIT('CC4HYPPRT called for an unknown Coupled '// 430 & 'Cluster model.') 431 END IF 432 433 IF (IPR4HYP.GT.5) THEN 434 WRITE(LUPRI,'(/1X,5(1X,A," operator",5X),4X,A,8X,A,/,121("-"))') 435 & 'A','B','C','D','E','delta','(asy. Resp.)' 436 ELSE 437 WRITE(LUPRI,'(/1X,5(1X,A," operator",5X),4X,A,/,86("-"))') 438 & 'A','B','C','D','E','delta' 439 END IF 440 441 ISYOLD = 1 442 443 DO IOPER = 1, N4ROPER 444 ISYTOT = 1 445 446 ISYMA = ISYOPR(IA4ROP(IOPER)) 447 ISYTOT = MULD2H(ISYTOT,ISYMA) 448 449 ISYMB = ISYOPR(IB4ROP(IOPER)) 450 ISYTOT = MULD2H(ISYTOT,ISYMB) 451 452 ISYMC = ISYOPR(IC4ROP(IOPER)) 453 ISYTOT = MULD2H(ISYTOT,ISYMC) 454 455 ISYMD = ISYOPR(ID4ROP(IOPER)) 456 ISYTOT = MULD2H(ISYTOT,ISYMD) 457 458 ISYME = ISYOPR(IE4ROP(IOPER)) 459 ISYTOT = MULD2H(ISYTOT,ISYME) 460 461 462 IFREQ = 1 463 IF (ISYTOT.EQ.1) THEN 464 IF (IPR4HYP.GT.5) THEN 465 WRITE(LUPRI,'(/1X,5(A7,F7.4,2X),G16.8," (",G10.3,")")') 466 & LBLOPR(IA4ROP(IOPER))(1:7),A4RFR(IFREQ), 467 & LBLOPR(IB4ROP(IOPER))(1:7),B4RFR(IFREQ), 468 & LBLOPR(IC4ROP(IOPER))(1:7),C4RFR(IFREQ), 469 & LBLOPR(ID4ROP(IOPER))(1:7),D4RFR(IFREQ), 470 & LBLOPR(IE4ROP(IOPER))(1:7),E4RFR(IFREQ), 471 & HALF*(HYPERPOL(IFREQ,IOPER,1)+HYPERPOL(IFREQ,IOPER,2)), 472 & HALF*(HYPERPOL(IFREQ,IOPER,1)-HYPERPOL(IFREQ,IOPER,2)) 473 ELSE 474 WRITE(LUPRI,'(/1X,5(A7,F7.4,2X),G16.8)') 475 & LBLOPR(IA4ROP(IOPER))(1:7),A4RFR(IFREQ), 476 & LBLOPR(IB4ROP(IOPER))(1:7),B4RFR(IFREQ), 477 & LBLOPR(IC4ROP(IOPER))(1:7),C4RFR(IFREQ), 478 & LBLOPR(ID4ROP(IOPER))(1:7),D4RFR(IFREQ), 479 & LBLOPR(IE4ROP(IOPER))(1:7),E4RFR(IFREQ), 480 & HALF*(HYPERPOL(IFREQ,IOPER,1)+HYPERPOL(IFREQ,IOPER,2)) 481 END IF 482 ISYOLD = 1 483 ELSE 484 IF (IPR4HYP.GT.5) THEN 485 IF (ISYOLD.EQ.1) WRITE(LUPRI,*) 486 WRITE(LUPRI,'(1X,5(A7,A7,2X),6X,A,7X," (",5X,A,6X,")")') 487 & LBLOPR(IA4ROP(IOPER))(1:7),' -.- ', 488 & LBLOPR(IB4ROP(IOPER))(1:7),' -.- ', 489 & LBLOPR(IC4ROP(IOPER))(1:7),' -.- ', 490 & LBLOPR(ID4ROP(IOPER))(1:7),' -.- ', 491 & LBLOPR(IE4ROP(IOPER))(1:7),' -.- ', 492 & '---', 493 & '---' 494 ELSE IF (IPR4HYP.GT.0) THEN 495 IF (ISYOLD.EQ.1) WRITE(LUPRI,*) 496 WRITE(LUPRI,'(1X,5(A7,A7,2X),6X,A,7X)') 497 & LBLOPR(IA4ROP(IOPER))(1:7),' -.- ', 498 & LBLOPR(IB4ROP(IOPER))(1:7),' -.- ', 499 & LBLOPR(IC4ROP(IOPER))(1:7),' -.- ', 500 & LBLOPR(ID4ROP(IOPER))(1:7),' -.- ', 501 & LBLOPR(IE4ROP(IOPER))(1:7),' -.- ', 502 & '---' 503 END IF 504 ISYOLD = 0 505 END IF 506 507 DO IFREQ = 2, N4RFREQ 508 IF (ISYTOT.EQ.1) THEN 509 IF (IPR4HYP.GT.5) THEN 510 WRITE(LUPRI,'(1X,5(7X,F7.4,2X),G16.8," (",G10.3,")")') 511 & A4RFR(IFREQ), B4RFR(IFREQ), C4RFR(IFREQ), D4RFR(IFREQ), 512 & E4RFR(IFREQ), 513 & HALF*(HYPERPOL(IFREQ,IOPER,1)+HYPERPOL(IFREQ,IOPER,2)), 514 & HALF*(HYPERPOL(IFREQ,IOPER,1)-HYPERPOL(IFREQ,IOPER,2)) 515 ELSE 516 WRITE(LUPRI,'(1X,5(7X,F7.4,2X),G16.8)') 517 & A4RFR(IFREQ), B4RFR(IFREQ), C4RFR(IFREQ), D4RFR(IFREQ), 518 & E4RFR(IFREQ), 519 & HALF*(HYPERPOL(IFREQ,IOPER,1)+HYPERPOL(IFREQ,IOPER,2)) 520 END IF 521 END IF 522 END DO 523 524 END DO 525 526 RETURN 527 END 528*---------------------------------------------------------------------* 529* END OF SUBROUTINE CC4HYPPRT * 530*---------------------------------------------------------------------* 531c /* deck CC4R_SETUP */ 532*=====================================================================* 533 SUBROUTINE CC4R_SETUP(MXTRAN3, MXVEC2, MXTRAN4, MXVEC1, 534 & I0HTRAN, I0HDOTS, W0H, N0HTRAN, 535 & I1HTRAN, I1HDOTS, W1H, N1HTRAN, 536 & I0GTRAN, I0GDOTS, W0G, N0GTRAN, 537 & I1GTRAN, I1GDOTS, W1G, N1GTRAN, 538 & I2GTRAN, I2GDOTS, W2G, N2GTRAN, 539 & I1FTRAN, I1FDOTS, W1F, N1FTRAN, 540 & I2FTRAN, I2FDOTS, W2F, N2FTRAN, 541 & I0FATRAN,I0FADOTS,W0FA,N0FATRAN, 542 & I1FATRAN,I1FADOTS,W1FA,N1FATRAN, 543 & I2FATRAN,I2FADOTS,W2FA,N2FATRAN, 544 & I2EATRAN,I2EADOTS,W2EA,N2EATRAN, 545 & IXTRAN, IXDOTS, WX, NXTRAN, 546 & HYPPOL, MXHYPPOL,LADD ) 547*---------------------------------------------------------------------* 548* 549* Purpose: set up for CC4R section 550* - list of H^0 matrix transformations 551* - list of H^A matrix transformations 552* - list of G^0 matrix transformations 553* - list of G^A matrix transformations 554* - list of G^AB matrix transformations 555* - list of F^A matrix transformations 556* - list of F^AB matrix transformations 557* - list of F^0{O} matrix transformations 558* - list of F^A{O} matrix transformations 559* - list of F^AB{O} matrix transformations 560* - list of ETA^AB{O} vector calculations 561* - list of chi vector dot products 562* 563* LADD = .FALSE. --> set lists 564* LADD = .TRUE. --> add contributions to hyperpolarizabilities 565* 566* Written by Christof Haettig, april 1997. 567* 568*=====================================================================* 569#if defined (IMPLICIT_NONE) 570 IMPLICIT NONE 571#else 572# include "implicit.h" 573#endif 574#include "priunit.h" 575#include "ccorb.h" 576#include "cc4rinf.h" 577#include "cc4perm.h" 578#include "ccroper.h" 579 580* local parameters: 581 CHARACTER*(20) MSGDBG 582 PARAMETER (MSGDBG = '[debug] CC4R_SETUP> ') 583 LOGICAL LOCDBG 584 PARAMETER (LOCDBG = .FALSE.) 585 586 LOGICAL LADD 587 588 INTEGER MXVEC2, MXTRAN3, MXVEC1, MXTRAN4, MXHYPPOL 589 590 INTEGER I0HTRAN(5,MXTRAN3) 591 INTEGER I0HDOTS(MXVEC2,MXTRAN3) 592 INTEGER I1HTRAN(5,MXTRAN4) 593 INTEGER I1HDOTS(MXVEC1,MXTRAN4) 594 595 INTEGER I0GTRAN(4,MXTRAN3) 596 INTEGER I0GDOTS(MXVEC2,MXTRAN3) 597 INTEGER I1GTRAN(4,MXTRAN3) 598 INTEGER I1GDOTS(MXVEC2,MXTRAN3) 599 INTEGER I2GTRAN(4,MXTRAN4) 600 INTEGER I2GDOTS(MXVEC1,MXTRAN4) 601 602 INTEGER I1FTRAN(3,MXTRAN3) 603 INTEGER I1FDOTS(MXVEC2,MXTRAN3) 604 INTEGER I2FTRAN(3,MXTRAN3) 605 INTEGER I2FDOTS(MXVEC2,MXTRAN3) 606 607 INTEGER I0FATRAN(5,MXTRAN3) 608 INTEGER I0FADOTS(MXVEC2,MXTRAN3) 609 INTEGER I1FATRAN(5,MXTRAN3) 610 INTEGER I1FADOTS(MXVEC2,MXTRAN3) 611 INTEGER I2FATRAN(5,MXTRAN4) 612 INTEGER I2FADOTS(MXVEC1,MXTRAN4) 613 614 INTEGER I2EATRAN(3,MXTRAN3) 615 INTEGER I2EADOTS(MXVEC2,MXTRAN3) 616 617 INTEGER IXTRAN(MXTRAN3) 618 INTEGER IXDOTS(MXVEC2,MXTRAN3) 619 620 INTEGER N0HTRAN, N0GTRAN, N0FATRAN 621 INTEGER N1HTRAN, N1GTRAN, N1FTRAN, N1FATRAN 622 INTEGER N2GTRAN, N2FTRAN, N2FATRAN, N2EATRAN, NXTRAN 623 INTEGER N4RRESF 624 INTEGER MXVH0, MXVH1, MXVG0, MXVG1, MXVG2, MXVF1, MXVF2 625 INTEGER MXVFA0, MXVFA1, MXVFA2, MXVEA2, MXX 626 627 INTEGER IVEC, ITRAN, I, IDX 628 INTEGER ISYML, ISYM1, ISYM2, ISYTOT, ISYM3 629 INTEGER IFREQ, IOPER 630 INTEGER P, ISIGN 631 632#if defined (SYS_CRAY) 633 REAL HYPPOL(MXHYPPOL) 634 REAL SIGN, SUM, SSUM 635 REAL H0CON(10), H1CON(5) 636 REAL G0CON(15), G1CON(30), G2CON(10) 637 REAL F1CON(15), F2CON(30) 638 REAL F0ACON(15), F1ACON(60), F2ACON(30) 639 REAL E2ACON(30) 640 REAL XCON(10) 641 REAL W0H(MXVEC2,MXTRAN3) 642 REAL W1H(MXVEC1,MXTRAN4) 643 REAL W0G(MXVEC2,MXTRAN3) 644 REAL W1G(MXVEC2,MXTRAN3) 645 REAL W2G(MXVEC1,MXTRAN4) 646 REAL W1F(MXVEC2,MXTRAN3) 647 REAL W2F(MXVEC2,MXTRAN3) 648 REAL W0FA(MXVEC2,MXTRAN3) 649 REAL W1FA(MXVEC2,MXTRAN3) 650 REAL W2FA(MXVEC1,MXTRAN4) 651 REAL W2EA(MXVEC2,MXTRAN3) 652 REAL WX(MXVEC2,MXTRAN3) 653#else 654 DOUBLE PRECISION SIGN, SUM, SSUM 655 DOUBLE PRECISION HYPPOL(MXHYPPOL) 656 DOUBLE PRECISION H0CON(10), H1CON(5) 657 DOUBLE PRECISION G0CON(15), G1CON(30), G2CON(10) 658 DOUBLE PRECISION F1CON(15), F2CON(30) 659 DOUBLE PRECISION F0ACON(15), F1ACON(60), F2ACON(30) 660 DOUBLE PRECISION E2ACON(30) 661 DOUBLE PRECISION XCON(10) 662 DOUBLE PRECISION W0H(MXVEC2,MXTRAN3) 663 DOUBLE PRECISION W1H(MXVEC1,MXTRAN4) 664 DOUBLE PRECISION W0G(MXVEC2,MXTRAN3) 665 DOUBLE PRECISION W1G(MXVEC2,MXTRAN3) 666 DOUBLE PRECISION W2G(MXVEC1,MXTRAN4) 667 DOUBLE PRECISION W1F(MXVEC2,MXTRAN3) 668 DOUBLE PRECISION W2F(MXVEC2,MXTRAN3) 669 DOUBLE PRECISION W0FA(MXVEC2,MXTRAN3) 670 DOUBLE PRECISION W1FA(MXVEC2,MXTRAN3) 671 DOUBLE PRECISION W2FA(MXVEC1,MXTRAN4) 672 DOUBLE PRECISION W2EA(MXVEC2,MXTRAN3) 673 DOUBLE PRECISION WX(MXVEC2,MXTRAN3) 674#endif 675 676 677 INTEGER IOP(5), ISY(5), IL1(5), IR1(5), IR2(10), IL2(10) 678 INTEGER IX3(10) 679 680 681* external functions: 682 INTEGER IR2TAMP 683 INTEGER IR1TAMP 684 INTEGER IL1ZETA 685 INTEGER IL2ZETA 686 INTEGER ICHI3 687 688 689*---------------------------------------------------------------------* 690* initializations: 691*---------------------------------------------------------------------* 692 IF (.NOT. LADD) THEN 693 N0HTRAN = 0 694 N1HTRAN = 0 695 N0GTRAN = 0 696 N1GTRAN = 0 697 N2GTRAN = 0 698 N1FTRAN = 0 699 N2FTRAN = 0 700 N0FATRAN = 0 701 N1FATRAN = 0 702 N2FATRAN = 0 703 N2EATRAN = 0 704 NXTRAN = 0 705 END IF 706 707 MXVH0 = 0 708 MXVH1 = 0 709 MXVG0 = 0 710 MXVG1 = 0 711 MXVG2 = 0 712 MXVF1 = 0 713 MXVF2 = 0 714 MXVFA0 = 0 715 MXVFA1 = 0 716 MXVFA2 = 0 717 MXVEA2 = 0 718 MXX = 0 719 720 N4RRESF = 0 721 722 CALL DZERO(HYPPOL,MXHYPPOL) 723 724*---------------------------------------------------------------------* 725* start loop over all requested third hyperpolarizabilities 726*---------------------------------------------------------------------* 727 728 DO IOPER = 1, N4ROPER 729 IOP(A) = IA4ROP(IOPER) 730 IOP(B) = IB4ROP(IOPER) 731 IOP(C) = IC4ROP(IOPER) 732 IOP(D) = ID4ROP(IOPER) 733 IOP(E) = IE4ROP(IOPER) 734 735 ISY(A) = ISYOPR(IOP(A)) 736 ISY(B) = ISYOPR(IOP(B)) 737 ISY(C) = ISYOPR(IOP(C)) 738 ISY(D) = ISYOPR(IOP(D)) 739 ISY(E) = ISYOPR(IOP(E)) 740 741 ISYTOT = 1 742 DO I = A, E 743 ISYTOT = MULD2H(ISYTOT,ISY(I)) 744 END DO 745 746 IF (ISYTOT .EQ. 1) THEN 747 748 DO IFREQ = 1, N4RFREQ 749 750 N4RRESF = N4RRESF + 1 751 752 DO ISIGN = 1, -1, -2 753 SIGN = DBLE(ISIGN) 754 755 IL1(A) =IL1ZETA(LBLOPR(IOP(A)),.FALSE.,SIGN*A4RFR(IFREQ),ISYML) 756 IL1(B) =IL1ZETA(LBLOPR(IOP(B)),.FALSE.,SIGN*B4RFR(IFREQ),ISYML) 757 IL1(C) =IL1ZETA(LBLOPR(IOP(C)),.FALSE.,SIGN*C4RFR(IFREQ),ISYML) 758 IL1(D) =IL1ZETA(LBLOPR(IOP(D)),.FALSE.,SIGN*D4RFR(IFREQ),ISYML) 759 IL1(E) =IL1ZETA(LBLOPR(IOP(E)),.FALSE.,SIGN*E4RFR(IFREQ),ISYML) 760 761 IR1(A) =IR1TAMP(LBLOPR(IOP(A)),.FALSE.,SIGN*A4RFR(IFREQ),ISYML) 762 IR1(B) =IR1TAMP(LBLOPR(IOP(B)),.FALSE.,SIGN*B4RFR(IFREQ),ISYML) 763 IR1(C) =IR1TAMP(LBLOPR(IOP(C)),.FALSE.,SIGN*C4RFR(IFREQ),ISYML) 764 IR1(D) =IR1TAMP(LBLOPR(IOP(D)),.FALSE.,SIGN*D4RFR(IFREQ),ISYML) 765 IR1(E) =IR1TAMP(LBLOPR(IOP(E)),.FALSE.,SIGN*E4RFR(IFREQ),ISYML) 766 767 IR2(AB)=IR2TAMP(LBLOPR(IOP(A)),.FALSE.,SIGN*A4RFR(IFREQ),ISYM1, 768 & LBLOPR(IOP(B)),.FALSE.,SIGN*B4RFR(IFREQ),ISYM2) 769 IR2(AC)=IR2TAMP(LBLOPR(IOP(A)),.FALSE.,SIGN*A4RFR(IFREQ),ISYM1, 770 & LBLOPR(IOP(C)),.FALSE.,SIGN*C4RFR(IFREQ),ISYM2) 771 IR2(AD)=IR2TAMP(LBLOPR(IOP(A)),.FALSE.,SIGN*A4RFR(IFREQ),ISYM1, 772 & LBLOPR(IOP(D)),.FALSE.,SIGN*D4RFR(IFREQ),ISYM2) 773 IR2(AE)=IR2TAMP(LBLOPR(IOP(A)),.FALSE.,SIGN*A4RFR(IFREQ),ISYM1, 774 & LBLOPR(IOP(E)),.FALSE.,SIGN*E4RFR(IFREQ),ISYM2) 775 IR2(BC)=IR2TAMP(LBLOPR(IOP(B)),.FALSE.,SIGN*B4RFR(IFREQ),ISYM1, 776 & LBLOPR(IOP(C)),.FALSE.,SIGN*C4RFR(IFREQ),ISYM2) 777 IR2(BD)=IR2TAMP(LBLOPR(IOP(B)),.FALSE.,SIGN*B4RFR(IFREQ),ISYM1, 778 & LBLOPR(IOP(D)),.FALSE.,SIGN*D4RFR(IFREQ),ISYM2) 779 IR2(BE)=IR2TAMP(LBLOPR(IOP(B)),.FALSE.,SIGN*B4RFR(IFREQ),ISYM1, 780 & LBLOPR(IOP(E)),.FALSE.,SIGN*E4RFR(IFREQ),ISYM2) 781 IR2(CD)=IR2TAMP(LBLOPR(IOP(C)),.FALSE.,SIGN*C4RFR(IFREQ),ISYM1, 782 & LBLOPR(IOP(D)),.FALSE.,SIGN*D4RFR(IFREQ),ISYM2) 783 IR2(CE)=IR2TAMP(LBLOPR(IOP(C)),.FALSE.,SIGN*C4RFR(IFREQ),ISYM1, 784 & LBLOPR(IOP(E)),.FALSE.,SIGN*E4RFR(IFREQ),ISYM2) 785 IR2(DE)=IR2TAMP(LBLOPR(IOP(D)),.FALSE.,SIGN*D4RFR(IFREQ),ISYM1, 786 & LBLOPR(IOP(E)),.FALSE.,SIGN*E4RFR(IFREQ),ISYM2) 787 788 IL2(AB) = IL2ZETA(LBLOPR(IOP(A)),SIGN*A4RFR(IFREQ),ISYM1, 789 & LBLOPR(IOP(B)),SIGN*B4RFR(IFREQ),ISYM2) 790 IL2(AC) = IL2ZETA(LBLOPR(IOP(A)),SIGN*A4RFR(IFREQ),ISYM1, 791 & LBLOPR(IOP(C)),SIGN*C4RFR(IFREQ),ISYM2) 792 IL2(AD) = IL2ZETA(LBLOPR(IOP(A)),SIGN*A4RFR(IFREQ),ISYM1, 793 & LBLOPR(IOP(D)),SIGN*D4RFR(IFREQ),ISYM2) 794 IL2(AE) = IL2ZETA(LBLOPR(IOP(A)),SIGN*A4RFR(IFREQ),ISYM1, 795 & LBLOPR(IOP(E)),SIGN*E4RFR(IFREQ),ISYM2) 796 IL2(BC) = IL2ZETA(LBLOPR(IOP(B)),SIGN*B4RFR(IFREQ),ISYM1, 797 & LBLOPR(IOP(C)),SIGN*C4RFR(IFREQ),ISYM2) 798 IL2(BD) = IL2ZETA(LBLOPR(IOP(B)),SIGN*B4RFR(IFREQ),ISYM1, 799 & LBLOPR(IOP(D)),SIGN*D4RFR(IFREQ),ISYM2) 800 IL2(BE) = IL2ZETA(LBLOPR(IOP(B)),SIGN*B4RFR(IFREQ),ISYM1, 801 & LBLOPR(IOP(E)),SIGN*E4RFR(IFREQ),ISYM2) 802 IL2(CD) = IL2ZETA(LBLOPR(IOP(C)),SIGN*C4RFR(IFREQ),ISYM1, 803 & LBLOPR(IOP(D)),SIGN*D4RFR(IFREQ),ISYM2) 804 IL2(CE) = IL2ZETA(LBLOPR(IOP(C)),SIGN*C4RFR(IFREQ),ISYM1, 805 & LBLOPR(IOP(E)),SIGN*E4RFR(IFREQ),ISYM2) 806 IL2(DE) = IL2ZETA(LBLOPR(IOP(D)),SIGN*D4RFR(IFREQ),ISYM1, 807 & LBLOPR(IOP(E)),SIGN*E4RFR(IFREQ),ISYM2) 808 809 IF (L_USE_CHI3) THEN 810 IX3(CDE) = ICHI3(LBLOPR(IOP(C)),SIGN*C4RFR(IFREQ),ISYM1, 811 & LBLOPR(IOP(D)),SIGN*D4RFR(IFREQ),ISYM2, 812 & LBLOPR(IOP(E)),SIGN*E4RFR(IFREQ),ISYM3) 813 IX3(BDE) = ICHI3(LBLOPR(IOP(B)),SIGN*B4RFR(IFREQ),ISYM1, 814 & LBLOPR(IOP(D)),SIGN*D4RFR(IFREQ),ISYM2, 815 & LBLOPR(IOP(E)),SIGN*E4RFR(IFREQ),ISYM3) 816 IX3(BCE) = ICHI3(LBLOPR(IOP(B)),SIGN*B4RFR(IFREQ),ISYM1, 817 & LBLOPR(IOP(C)),SIGN*C4RFR(IFREQ),ISYM2, 818 & LBLOPR(IOP(E)),SIGN*E4RFR(IFREQ),ISYM3) 819 IX3(BCD) = ICHI3(LBLOPR(IOP(B)),SIGN*B4RFR(IFREQ),ISYM1, 820 & LBLOPR(IOP(C)),SIGN*C4RFR(IFREQ),ISYM2, 821 & LBLOPR(IOP(D)),SIGN*D4RFR(IFREQ),ISYM3) 822 IX3(ADE) = ICHI3(LBLOPR(IOP(A)),SIGN*A4RFR(IFREQ),ISYM1, 823 & LBLOPR(IOP(D)),SIGN*D4RFR(IFREQ),ISYM2, 824 & LBLOPR(IOP(E)),SIGN*E4RFR(IFREQ),ISYM3) 825 IX3(ACE) = ICHI3(LBLOPR(IOP(A)),SIGN*A4RFR(IFREQ),ISYM1, 826 & LBLOPR(IOP(C)),SIGN*C4RFR(IFREQ),ISYM2, 827 & LBLOPR(IOP(E)),SIGN*E4RFR(IFREQ),ISYM3) 828 IX3(ACD) = ICHI3(LBLOPR(IOP(A)),SIGN*A4RFR(IFREQ),ISYM1, 829 & LBLOPR(IOP(C)),SIGN*C4RFR(IFREQ),ISYM2, 830 & LBLOPR(IOP(D)),SIGN*D4RFR(IFREQ),ISYM3) 831 IX3(ABE) = ICHI3(LBLOPR(IOP(A)),SIGN*A4RFR(IFREQ),ISYM1, 832 & LBLOPR(IOP(B)),SIGN*B4RFR(IFREQ),ISYM2, 833 & LBLOPR(IOP(E)),SIGN*E4RFR(IFREQ),ISYM3) 834 IX3(ABD) = ICHI3(LBLOPR(IOP(A)),SIGN*A4RFR(IFREQ),ISYM1, 835 & LBLOPR(IOP(B)),SIGN*B4RFR(IFREQ),ISYM2, 836 & LBLOPR(IOP(D)),SIGN*D4RFR(IFREQ),ISYM3) 837 IX3(ABC) = ICHI3(LBLOPR(IOP(A)),SIGN*A4RFR(IFREQ),ISYM1, 838 & LBLOPR(IOP(B)),SIGN*B4RFR(IFREQ),ISYM2, 839 & LBLOPR(IOP(C)),SIGN*C4RFR(IFREQ),ISYM3) 840 END IF 841*---------------------------------------------------------------------* 842* set up list of H^0 matrix transformations, 10 permutations 843*---------------------------------------------------------------------* 844 DO P = 1, 5 845 CALL CC_SETH1112(I0HTRAN,I0HDOTS,MXTRAN3,MXVEC2, 846 & 0,IR1(J3(P)),IR1(J4(P)),IR1(J5(P)),IR2(JP1(P)),ITRAN,IVEC) 847 N0HTRAN = MAX(N0HTRAN,ITRAN) 848 MXVH0 = MAX(MXVH0,IVEC) 849 H0CON(P) = W0H(IVEC,ITRAN) 850 851 CALL CC_SETH1112(I0HTRAN,I0HDOTS,MXTRAN3,MXVEC2, 852 & 0,IR1(J2(P)),IR1(J4(P)),IR1(J5(P)),IR2(JP2(P)),ITRAN,IVEC) 853 N0HTRAN = MAX(N0HTRAN,ITRAN) 854 MXVH0 = MAX(MXVH0,IVEC) 855 H0CON(P+5) = W0H(IVEC,ITRAN) 856 END DO 857 858*---------------------------------------------------------------------* 859* set up list of H^A matrix transformations, 5 permutations 860*---------------------------------------------------------------------* 861 DO P = 1, 5 862 CALL CC_SETH1111(I1HTRAN,I1HDOTS,MXTRAN4,MXVEC1, 863 & IL1(J1(P)),IR1(J2(P)),IR1(J3(P)),IR1(J4(P)),IR1(J5(P)), 864 & ITRAN,IVEC) 865 N1HTRAN = MAX(N1HTRAN,ITRAN) 866 MXVH1 = MAX(MXVH1,IVEC) 867 H1CON(P) = W1H(IVEC,ITRAN) 868 END DO 869 870*---------------------------------------------------------------------* 871* set up list of G^0 matrix transformations, 15 permutations 872*---------------------------------------------------------------------* 873 DO P = 1, 15 874 CALL CC_SETG212(I0GTRAN,I0GDOTS,MXTRAN3,MXVEC2, 875 & 0,IR2(IP1(P)),IR1(I1(P)),IR2(IP2(P)),ITRAN,IVEC) 876 N0GTRAN = MAX(N0GTRAN,ITRAN) 877 MXVG0 = MAX(MXVG0,IVEC) 878 G0CON(P) = W0G(IVEC,ITRAN) 879 END DO 880 881*---------------------------------------------------------------------* 882* set up list of G^A matrix transformations, 30 permutations 883*---------------------------------------------------------------------* 884 DO P = 1, 15 885 CALL CC_SETG112(I1GTRAN,I1GDOTS,MXTRAN3,MXVEC2, 886 & IL1(I1(P)),IR1(I2(P)),IR1(I3(P)),IR2(IP2(P)),ITRAN,IVEC) 887 N1GTRAN = MAX(N1GTRAN,ITRAN) 888 MXVG1 = MAX(MXVG1,IVEC) 889 G1CON(P) = W1G(IVEC,ITRAN) 890 891 CALL CC_SETG112(I1GTRAN,I1GDOTS,MXTRAN3,MXVEC2, 892 & IL1(I1(P)),IR1(I4(P)),IR1(I5(P)),IR2(IP1(P)),ITRAN,IVEC) 893 N1GTRAN = MAX(N1GTRAN,ITRAN) 894 MXVG1 = MAX(MXVG1,IVEC) 895 G1CON(P+15) = W1G(IVEC,ITRAN) 896 END DO 897 898*---------------------------------------------------------------------* 899* set up list of G^AB matrix transformations, 10 permutations 900*---------------------------------------------------------------------* 901 DO P = 1, 5 902 CALL CCQR_SETG(I2GTRAN,I2GDOTS,MXTRAN4,MXVEC1, 903 & IL2(JP1(P)),IR1(J3(P)),IR1(J4(P)),IR1(J5(P)),ITRAN,IVEC) 904 N2GTRAN = MAX(N2GTRAN,ITRAN) 905 MXVG2 = MAX(MXVG2,IVEC) 906 G2CON(P) = W2G(IVEC,ITRAN) 907 908 CALL CCQR_SETG(I2GTRAN,I2GDOTS,MXTRAN4,MXVEC1, 909 & IL2(JP2(P)),IR1(J2(P)),IR1(J4(P)),IR1(J5(P)),ITRAN,IVEC) 910 N2GTRAN = MAX(N2GTRAN,ITRAN) 911 MXVG2 = MAX(MXVG2,IVEC) 912 G2CON(P+5) = W2G(IVEC,ITRAN) 913 END DO 914 915*---------------------------------------------------------------------* 916* set up list of F^A matrix transformations, 15 permutations 917*---------------------------------------------------------------------* 918 DO P = 1, 15 919 CALL CCQR_SETF(I1FTRAN,I1FDOTS,MXTRAN3,MXVEC2, 920 & IL1(I1(P)),IR2(IP1(P)),IR2(IP2(P)),ITRAN,IVEC) 921 N1FTRAN = MAX(N1FTRAN,ITRAN) 922 MXVF1 = MAX(MXVF1,IVEC) 923 F1CON(P) = W1F(IVEC,ITRAN) 924 END DO 925 926*---------------------------------------------------------------------* 927* set up list of F^AB matrix transformations, 30 permutations 928*---------------------------------------------------------------------* 929 DO P = 1, 15 930 CALL CC_SETF12(I2FTRAN,I2FDOTS,MXTRAN3,MXVEC2, 931 & IL2(IP1(P)),IR1(I1(P)),IR2(IP2(P)),ITRAN,IVEC) 932 N2FTRAN = MAX(N2FTRAN,ITRAN) 933 MXVF2 = MAX(MXVF2,IVEC) 934 F2CON(P) = W2F(IVEC,ITRAN) 935 936 CALL CC_SETF12(I2FTRAN,I2FDOTS,MXTRAN3,MXVEC2, 937 & IL2(IP2(P)),IR1(I1(P)),IR2(IP1(P)),ITRAN,IVEC) 938 N2FTRAN = MAX(N2FTRAN,ITRAN) 939 MXVF2 = MAX(MXVF2,IVEC) 940 F2CON(P+15) = W2F(IVEC,ITRAN) 941 END DO 942 943*---------------------------------------------------------------------* 944* set up list of F^0{O} matrix transformations, 15 permutations 945*---------------------------------------------------------------------* 946 DO P = 1, 15 947 CALL CCQR_SETFA(I0FATRAN,I0FADOTS,MXTRAN3,MXVEC2, 948 & 0,IOP(I1(P)),IR2(IP1(P)),IR2(IP2(P)),ITRAN,IVEC) 949 N0FATRAN = MAX(N0FATRAN,ITRAN) 950 MXVFA0 = MAX(MXVFA0,IVEC) 951 F0ACON(P) = W0FA(IVEC,ITRAN) 952 END DO 953 954*---------------------------------------------------------------------* 955* set up list of F^A{O} matrix transformations, 60 permutations 956*---------------------------------------------------------------------* 957 DO P = 1, 15 958 CALL CC_SETFA12(I1FATRAN,I1FADOTS,MXTRAN3,MXVEC2, 959 & IL1(I1(P)),IOP(I2(P)),IR1(I3(P)),IR2(IP2(P)),ITRAN,IVEC) 960 N1FATRAN = MAX(N1FATRAN,ITRAN) 961 MXVFA1 = MAX(MXVFA1,IVEC) 962 F1ACON(P) = W1FA(IVEC,ITRAN) 963 964 CALL CC_SETFA12(I1FATRAN,I1FADOTS,MXTRAN3,MXVEC2, 965 & IL1(I1(P)),IOP(I3(P)),IR1(I2(P)),IR2(IP2(P)),ITRAN,IVEC) 966 N1FATRAN = MAX(N1FATRAN,ITRAN) 967 MXVFA1 = MAX(MXVFA1,IVEC) 968 F1ACON(P+15) = W1FA(IVEC,ITRAN) 969 970 CALL CC_SETFA12(I1FATRAN,I1FADOTS,MXTRAN3,MXVEC2, 971 & IL1(I1(P)),IOP(I4(P)),IR1(I5(P)),IR2(IP1(P)),ITRAN,IVEC) 972 N1FATRAN = MAX(N1FATRAN,ITRAN) 973 MXVFA1 = MAX(MXVFA1,IVEC) 974 F1ACON(P+30) = W1FA(IVEC,ITRAN) 975 976 CALL CC_SETFA12(I1FATRAN,I1FADOTS,MXTRAN3,MXVEC2, 977 & IL1(I1(P)),IOP(I5(P)),IR1(I4(P)),IR2(IP1(P)),ITRAN,IVEC) 978 N1FATRAN = MAX(N1FATRAN,ITRAN) 979 MXVFA1 = MAX(MXVFA1,IVEC) 980 F1ACON(P+45) = W1FA(IVEC,ITRAN) 981 END DO 982 983*---------------------------------------------------------------------* 984* set up list of F^AB{O} matrix transformations, 30 permutations 985*---------------------------------------------------------------------* 986 DO P = 1, 15 987 CALL CCQR_SETFA(I2FATRAN,I2FADOTS,MXTRAN4,MXVEC1, 988 & IL2(IP1(P)),IOP(I1(P)),IR1(I4(P)),IR1(I5(P)),ITRAN,IVEC) 989 N2FATRAN = MAX(N2FATRAN,ITRAN) 990 MXVFA2 = MAX(MXVFA2,IVEC) 991 F2ACON(P) = W2FA(IVEC,ITRAN) 992 993 CALL CCQR_SETFA(I2FATRAN,I2FADOTS,MXTRAN4,MXVEC1, 994 & IL2(IP2(P)),IOP(I1(P)),IR1(I2(P)),IR1(I3(P)),ITRAN,IVEC) 995 N2FATRAN = MAX(N2FATRAN,ITRAN) 996 MXVFA2 = MAX(MXVFA2,IVEC) 997 F2ACON(P+15) = W2FA(IVEC,ITRAN) 998 END DO 999 1000*---------------------------------------------------------------------* 1001* set up list of ETA{O} vector calculations, 30 permutations 1002*---------------------------------------------------------------------* 1003 DO P = 1, 15 1004 CALL CCQR_SETEA(I2EATRAN,I2EADOTS,MXTRAN3,MXVEC2, 1005 & IL2(IP1(P)),IOP(I1(P)),IR2(IP2(P)),ITRAN,IVEC) 1006 N2EATRAN = MAX(N2EATRAN,ITRAN) 1007 MXVEA2 = MAX(MXVEA2,IVEC) 1008 E2ACON(P) = W2EA(IVEC,ITRAN) 1009 1010 CALL CCQR_SETEA(I2EATRAN,I2EADOTS,MXTRAN3,MXVEC2, 1011 & IL2(IP2(P)),IOP(I1(P)),IR2(IP1(P)),ITRAN,IVEC) 1012 N2EATRAN = MAX(N2EATRAN,ITRAN) 1013 MXVEA2 = MAX(MXVEA2,IVEC) 1014 E2ACON(P+15) = W2EA(IVEC,ITRAN) 1015 END DO 1016 1017*---------------------------------------------------------------------* 1018* set up list of chi vector dot products, 6 permutations 1019*---------------------------------------------------------------------* 1020 IF (L_USE_CHI3) THEN 1021 DO P = 1, 10 1022 CALL CC_SETDOT(IXTRAN,IXDOTS,MXTRAN3,MXVEC2, 1023 & IX3(11-P),IR2(P),ITRAN,IVEC) 1024 NXTRAN = MAX(NXTRAN,ITRAN) 1025 MXX = MAX(MXX,IVEC) 1026 XCON(P) = WX(IVEC,ITRAN) 1027 END DO 1028 ELSE 1029 DO P = 1, 10 1030 XCON(P) = 0.0d0 1031 END DO 1032 END IF 1033 1034 1035*---------------------------------------------------------------------* 1036* add all 250 contributions to the hyperpolarizabilities: 1037*---------------------------------------------------------------------* 1038 IF (LADD) THEN 1039 IDX = N4RFREQ*(IOPER-1) + IFREQ 1040 IF (ISIGN.EQ.-1) IDX = IDX + N4RFREQ*N4ROPER 1041 1042 DO P = 1, 5 1043 HYPPOL(IDX) = HYPPOL(IDX) + H0CON(P) + H0CON(P+5) + H1CON(P) + 1044 & G2CON(P) + G2CON(P+5) + 1045 & XCON(P) + XCON(P+5) 1046 END DO 1047 1048 DO P = 1, 15 1049 HYPPOL(IDX) = HYPPOL(IDX) + 1050 & G0CON(P) + G1CON(P) + G1CON(P+15) + 1051 & F1CON(P) + F2CON(P) + F2CON(P+15) + F0ACON(P) + 1052 & F1ACON(P) + F1ACON(P+15) + F1ACON(P+30) + F1ACON(P+45) + 1053 & F2ACON(P) + F2ACON(P+15) + E2ACON(P) + E2ACON(P+15) 1054 END DO 1055 1056 IF (LOCDBG) THEN 1057 WRITE (LUPRI,*) 'IOP:',IOP 1058 WRITE (LUPRI,*) 'IDX:',IDX 1059 1060 SSUM= 0.0d0 1061 SUM = 0.0d0 1062 DO I = 1, 10 1063 SUM = SUM + H0CON(I) 1064 END DO 1065 SSUM = SSUM + SUM 1066 WRITE (LUPRI,*) 'H0CON:',H0CON 1067 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1068 1069 SUM = 0.0d0 1070 DO I = 1, 5 1071 SUM = SUM + H1CON(I) 1072 END DO 1073 SSUM = SSUM + SUM 1074 WRITE (LUPRI,*) 'H1CON:',H1CON 1075 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1076 1077 1078 SUM = 0.0d0 1079 DO I = 1, 15 1080 SUM = SUM + G0CON(I) 1081 END DO 1082 SSUM = SSUM + SUM 1083 WRITE (LUPRI,*) 'G0CON:',G0CON 1084 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1085 1086 1087 SUM = 0.0d0 1088 DO I = 1, 30 1089 SUM = SUM + G1CON(I) 1090 END DO 1091 SSUM = SSUM + SUM 1092 WRITE (LUPRI,*) 'G1CON:',G1CON 1093 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1094 1095 SUM = 0.0d0 1096 DO I = 1, 10 1097 SUM = SUM + G2CON(I) 1098 END DO 1099 SSUM = SSUM + SUM 1100 WRITE (LUPRI,*) 'G2CON:',G2CON 1101 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1102 1103 SUM = 0.0d0 1104 DO I = 1, 15 1105 SUM = SUM + F1CON(I) 1106 END DO 1107 SSUM = SSUM + SUM 1108 WRITE (LUPRI,*) 'F1CON:',F1CON 1109 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1110 1111 SUM = 0.0d0 1112 DO I = 1, 30 1113 SUM = SUM + F2CON(I) 1114 END DO 1115 SSUM = SSUM + SUM 1116 WRITE (LUPRI,*) 'F2CON:',F2CON 1117 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1118 1119 1120 SUM = 0.0d0 1121 DO I = 1, 15 1122 SUM = SUM + F0ACON(I) 1123 END DO 1124 SSUM = SSUM + SUM 1125 WRITE (LUPRI,*) 'F0ACON:',F0ACON 1126 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1127 1128 1129 SUM = 0.0d0 1130 DO I = 1, 60 1131 SUM = SUM + F1ACON(I) 1132 END DO 1133 SSUM = SSUM + SUM 1134 WRITE (LUPRI,*) 'F1ACON:',F1ACON 1135 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1136 1137 1138 SUM = 0.0d0 1139 DO I = 1, 30 1140 SUM = SUM + F2ACON(I) 1141 END DO 1142 SSUM = SSUM + SUM 1143 WRITE (LUPRI,*) 'F2ACON:',F2ACON 1144 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1145 1146 1147 SUM = 0.0d0 1148 DO I = 1, 30 1149 SUM = SUM + E2ACON(I) 1150 END DO 1151 SSUM = SSUM + SUM 1152 WRITE (LUPRI,*) 'E2ACON:',E2ACON 1153 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1154 1155 SUM = 0.0d0 1156 DO I = 1, 30 1157 SUM = SUM + XCON(I) 1158 END DO 1159 SSUM = SSUM + SUM 1160 WRITE (LUPRI,*) 'XCON:',XCON 1161 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1162 1163 WRITE (LUPRI,*) 'HYPPOL:',HYPPOL(IDX) 1164 END IF 1165 1166 END IF 1167 1168*---------------------------------------------------------------------* 1169* end loop over all requested hyperpolarizabilities 1170*---------------------------------------------------------------------* 1171 END DO 1172 END DO 1173 END IF 1174 END DO 1175 1176*---------------------------------------------------------------------* 1177* print statistics and lists: 1178*---------------------------------------------------------------------* 1179 IF (.NOT. LADD) THEN 1180 1181* general statistics: 1182 WRITE(LUPRI,'(/,/3X,A,I4,A)') 'For the requested',N4RRESF, 1183 & ' quartic response functions ' 1184 WRITE(LUPRI,'((8X,A,I3,A))') 1185 & ' - ',N0HTRAN, ' H^0 matrix transformations ', 1186 & ' - ',N1HTRAN, ' H^A matrix transformations ', 1187 & ' - ',N0GTRAN, ' G^0 matrix transformations ', 1188 & ' - ',N1GTRAN, ' G^A matrix transformations ', 1189 & ' - ',N2GTRAN, ' G^AB matrix transformations ', 1190 & ' - ',N1FTRAN, ' F^A matrix transformations ', 1191 & ' - ',N2FTRAN, ' F^AB matrix transformations ', 1192 & ' - ',N0FATRAN, ' F^0{O} matrix transformations ', 1193 & ' - ',N1FATRAN, ' F^A{O} matrix transformations ', 1194 & ' - ',N2FATRAN, ' F^AB{O} matrix transformations ', 1195 & ' - ',N2EATRAN, ' ETA^AB{O} vector calculations ' 1196 WRITE(LUPRI,'(3X,A,/,/)') 'will be performed.' 1197 1198 1199* H^0 matrix transformations: 1200 IF (LOCDBG) THEN 1201 WRITE (LUPRI,*) 'List of H^0 matrix transformations:' 1202 DO ITRAN = 1, N0HTRAN 1203 WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG, 1204 & (I0HTRAN(I,ITRAN),I=1,4),(I0HDOTS(I,ITRAN),I=1,MXVH0) 1205 END DO 1206 WRITE (LUPRI,*) 1207 END IF 1208 1209* H^A matrix transformations: 1210 IF (LOCDBG) THEN 1211 WRITE (LUPRI,*) 'List of H^A matrix transformations:' 1212 DO ITRAN = 1, N1HTRAN 1213 WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG, 1214 & (I1HTRAN(I,ITRAN),I=1,4),(I1HDOTS(I,ITRAN),I=1,MXVH1) 1215 END DO 1216 WRITE (LUPRI,*) 1217 END IF 1218 1219* G^0 matrix transformations: 1220 IF (LOCDBG) THEN 1221 WRITE (LUPRI,*) 'List of G^0 matrix transformations:' 1222 DO ITRAN = 1, N0GTRAN 1223 WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG, 1224 & (I0GTRAN(I,ITRAN),I=1,3),(I0GDOTS(I,ITRAN),I=1,MXVG0) 1225 END DO 1226 WRITE (LUPRI,*) 1227 END IF 1228 1229* G^A matrix transformations: 1230 IF (LOCDBG) THEN 1231 WRITE (LUPRI,*) 'List of G^A matrix transformations:' 1232 DO ITRAN = 1, N1GTRAN 1233 WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG, 1234 & (I1GTRAN(I,ITRAN),I=1,3),(I1GDOTS(I,ITRAN),I=1,MXVG1) 1235 END DO 1236 WRITE (LUPRI,*) 1237 END IF 1238 1239* G^AB matrix transformations: 1240 IF (LOCDBG) THEN 1241 WRITE (LUPRI,*) 'List of G^AB matrix transformations:' 1242 DO ITRAN = 1, N2GTRAN 1243 WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG, 1244 & (I2GTRAN(I,ITRAN),I=1,3),(I2GDOTS(I,ITRAN),I=1,MXVG2) 1245 END DO 1246 WRITE (LUPRI,*) 1247 END IF 1248 1249* F^A matrix transformations: 1250 IF (LOCDBG) THEN 1251 WRITE (LUPRI,*) 'List of F^A matrix transformations:' 1252 DO ITRAN = 1, N1FTRAN 1253 WRITE(LUPRI,'(A,2I5,5X,(12I5,20X))') MSGDBG, 1254 & (I1FTRAN(I,ITRAN),I=1,2),(I1FDOTS(I,ITRAN),I=1,MXVF1) 1255 END DO 1256 WRITE (LUPRI,*) 1257 END IF 1258 1259* F^AB matrix transformations: 1260 IF (LOCDBG) THEN 1261 WRITE (LUPRI,*) 'List of F^AB matrix transformations:' 1262 DO ITRAN = 1, N2FTRAN 1263 WRITE(LUPRI,'(A,2I5,5X,(12I5,20X))') MSGDBG, 1264 & (I2FTRAN(I,ITRAN),I=1,2),(I2FDOTS(I,ITRAN),I=1,MXVF2) 1265 END DO 1266 WRITE (LUPRI,*) 1267 END IF 1268 1269* F^0{O} matrix transformations: 1270 IF (LOCDBG) THEN 1271 WRITE (LUPRI,*) 'List of F{O} matrix transformations:' 1272 DO ITRAN = 1, N0FATRAN 1273 WRITE(LUPRI,'(A,5I5,5X,(12I5,20X))') MSGDBG, 1274 & (I0FATRAN(I,ITRAN),I=1,5),(I0FADOTS(I,ITRAN),I=1,MXVFA0) 1275 END DO 1276 WRITE (LUPRI,*) 1277 END IF 1278 1279* F^A{O} matrix transformations: 1280 IF (LOCDBG) THEN 1281 WRITE (LUPRI,*) 'List of F^A{O} matrix transformations:' 1282 DO ITRAN = 1, N1FATRAN 1283 WRITE(LUPRI,'(A,5I5,5X,(12I5,20X))') MSGDBG, 1284 & (I1FATRAN(I,ITRAN),I=1,5),(I1FADOTS(I,ITRAN),I=1,MXVFA1) 1285 END DO 1286 WRITE (LUPRI,*) 1287 END IF 1288 1289* F^AB{O} matrix transformations: 1290 IF (LOCDBG) THEN 1291 WRITE (LUPRI,*) 'List of F^AB{O} matrix transformations:' 1292 DO ITRAN = 1, N2FATRAN 1293 WRITE(LUPRI,'(A,5I5,5X,(12I5,20X))') MSGDBG, 1294 & (I2FATRAN(I,ITRAN),I=1,5),(I2FADOTS(I,ITRAN),I=1,MXVFA2) 1295 END DO 1296 WRITE (LUPRI,*) 1297 END IF 1298 1299* ETA^AB{O} vector calculations: 1300 IF (LOCDBG) THEN 1301 WRITE (LUPRI,*) 'List of ETA^AB{O} vector calculations:' 1302 DO ITRAN = 1, N2EATRAN 1303 WRITE(LUPRI,'(A,2I5,5X,(12I5,20X))') MSGDBG, 1304 & (I2EATRAN(I,ITRAN),I=1,2),(I2EADOTS(I,ITRAN),I=1,MXVEA2) 1305 END DO 1306 WRITE (LUPRI,*) 1307 END IF 1308 1309 CALL FLSHFO(LUPRI) 1310 1311 END IF ! (.NOT. LADD) 1312*---------------------------------------------------------------------* 1313* return 1314*---------------------------------------------------------------------* 1315 1316 RETURN 1317 END 1318 1319*---------------------------------------------------------------------* 1320* END OF SUBROUTINE CC4R_SETUP * 1321*---------------------------------------------------------------------* 1322 1323c /* deck CC_SETH1112 */ 1324*=====================================================================* 1325 SUBROUTINE CC_SETH1112(IHTRAN,IHDOTS,MXTRAN,MXVEC,IZETAV, 1326 & ITAMPA,ITAMPB,ITAMPC,ITAMPD,ITRAN,IVEC) 1327*---------------------------------------------------------------------* 1328* 1329* Purpose: maintain a list of dot products of H matrix 1330* transformations with 3 right amplitude vectors: 1331* (Z*K*T^A*T^B*T^C) * T^D 1332* assumes that T^A, T^B and T^C belong to the same list 1333* and T^D belongs to a second list 1334* 1335* IHTRAN - list of H matrix transformations 1336* IHDOTS - list of vectors it should be dottet on 1337* 1338* MXTRAN - maximum list dimension 1339* MXVEC - maximum second dimesion for IHDOTS 1340* 1341* IZETAV - index of lagrangian multiplier vector 1342* ITAMPA - index of amplitude vector A 1343* ITAMPB - index of amplitude vector B 1344* ITAMPC - index of amplitude vector C 1345* ITAMPD - index of amplitude vector D 1346* 1347* ITRAN - index in IHTRAN list 1348* IVEC - second index in IHDOTS list 1349* 1350* Written by Christof Haettig, april 1997. 1351* 1352*=====================================================================* 1353 IMPLICIT NONE 1354 INTEGER MXVEC, MXTRAN 1355 INTEGER IHTRAN(5,MXTRAN) 1356 INTEGER IHDOTS(MXVEC,MXTRAN) 1357 1358 LOGICAL LFNDA, LFNDB, LFNDC, LFNDD 1359 INTEGER IZETAV, ITAMPA, ITAMPB, ITAMPC, ITAMPD 1360 INTEGER ITRAN, IVEC 1361 INTEGER ITAMP, I, IDX 1362 1363* statement functions: 1364#include "priunit.h" 1365 LOGICAL LHTST, LHEND 1366 INTEGER IL, IA, IB,IC 1367 LHTST(ITRAN,IL,IA,IB,IC) = IHTRAN(1,ITRAN).EQ.IL .AND. ( 1368 & (IHTRAN(2,ITRAN).EQ.IA .AND. IHTRAN(3,ITRAN).EQ.IB 1369 & .AND. IHTRAN(4,ITRAN).EQ.IC) .OR. 1370 & (IHTRAN(2,ITRAN).EQ.IB .AND. IHTRAN(3,ITRAN).EQ.IA 1371 & .AND. IHTRAN(4,ITRAN).EQ.IC) .OR. 1372 & (IHTRAN(2,ITRAN).EQ.IC .AND. IHTRAN(3,ITRAN).EQ.IA 1373 & .AND. IHTRAN(4,ITRAN).EQ.IB) .OR. 1374 & (IHTRAN(2,ITRAN).EQ.IA .AND. IHTRAN(3,ITRAN).EQ.IC 1375 & .AND. IHTRAN(4,ITRAN).EQ.IB) .OR. 1376 & (IHTRAN(2,ITRAN).EQ.IB .AND. IHTRAN(3,ITRAN).EQ.IC 1377 & .AND. IHTRAN(4,ITRAN).EQ.IA) .OR. 1378 & (IHTRAN(2,ITRAN).EQ.IC .AND. IHTRAN(3,ITRAN).EQ.IB 1379 & .AND. IHTRAN(4,ITRAN).EQ.IA) ) 1380 LHEND(ITRAN) = ITRAN.GT.MXTRAN .OR. 1381 & (IHTRAN(1,ITRAN)+IHTRAN(2,ITRAN)+ 1382 & IHTRAN(3,ITRAN)+IHTRAN(4,ITRAN) ).LE.0 1383 1384*---------------------------------------------------------------------* 1385* set up list of K matrix transformations 1386*---------------------------------------------------------------------* 1387 ITRAN = 1 1388 LFNDD = LHTST(ITRAN,IZETAV,ITAMPA,ITAMPB,ITAMPC) 1389 1390 DO WHILE ( .NOT. (LFNDD.OR.LHEND(ITRAN)) ) 1391 ITRAN = ITRAN + 1 1392 LFNDD = LHTST(ITRAN,IZETAV,ITAMPA,ITAMPB,ITAMPC) 1393 END DO 1394 1395 IF (.NOT.LFNDD) THEN 1396 IHTRAN(1,ITRAN) = IZETAV 1397 IHTRAN(2,ITRAN) = ITAMPA 1398 IHTRAN(3,ITRAN) = ITAMPB 1399 IHTRAN(4,ITRAN) = ITAMPC 1400 ITAMP = ITAMPD 1401 ELSE 1402 ITAMP = ITAMPD 1403 END IF 1404 1405 IVEC = 1 1406 DO WHILE (IHDOTS(IVEC,ITRAN).NE.ITAMP .AND. 1407 & IHDOTS(IVEC,ITRAN).NE.0 .AND. IVEC.LE.MXVEC) 1408 IVEC = IVEC + 1 1409 END DO 1410 1411 IHDOTS(IVEC,ITRAN) = ITAMP 1412 1413*---------------------------------------------------------------------* 1414 IF (IVEC.GT.MXVEC .OR. ITRAN.GT.MXTRAN) THEN 1415 WRITE (LUPRI,*) 'IVEC :',IVEC 1416 WRITE (LUPRI,*) 'ITRAN:',ITRAN 1417 WRITE (LUPRI,*) 'ITAMPA,ITAMPB,ITAMPC,ITAMPD:', 1418 & ITAMPA,ITAMPB,ITAMPC,ITAMPD 1419 IDX = 1 1420 DO WHILE (.NOT. LHEND(IDX)) 1421 WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') 'CC_SETH1112>', 1422 & (IHTRAN(I,IDX),I=1,4),(IHDOTS(I,IDX),I=1,MXVEC) 1423 IDX = IDX + 1 1424 END DO 1425 CALL FLSHFO(LUPRI) 1426 CALL QUIT('Overflow error in CC_SETH1112.') 1427 END IF 1428 1429 RETURN 1430 END 1431*---------------------------------------------------------------------* 1432* END OF SUBROUTINE CC_SETH1112 * 1433*---------------------------------------------------------------------* 1434c /* deck CC_SETG212 */ 1435*=====================================================================* 1436 SUBROUTINE CC_SETG212(IGTRAN,IGDOTS,MXTRAN,MXVEC, 1437 & IZETAV,ITAMPA,ITAMPB,ITAMPC,ITRAN,IVEC) 1438*---------------------------------------------------------------------* 1439* 1440* Purpose: maintain a list of dot products of G matrix 1441* transformations with right amplitude vectors: 1442* (Z*G*T^A*T^B) * T^C 1443* assumes that T^A and T^C belong to one list, 1444* and T^B belongs to a second list 1445* 1446* IGTRAN - list of G matrix transformations 1447* IGDOTS - list of vectors it should be dottet on 1448* 1449* MXTRAN - maximum list dimension 1450* MXVEC - maximum second dimesion for IGDOTS 1451* 1452* IZETAV - index of lagrangian multiplier vector 1453* ITAMPA - index of amplitude vector A 1454* ITAMPB - index of amplitude vector B 1455* ITAMPC - index of amplitude vector C 1456* 1457* ITRAN - index in IGTRAN list 1458* IVEC - second index in IGDOTS list 1459* 1460* Written by Christof Haettig, april 1997. 1461* 1462*=====================================================================* 1463 IMPLICIT NONE 1464 1465 LOGICAL LOCDBG 1466 PARAMETER (LOCDBG = .FALSE.) 1467 1468 INTEGER MXVEC, MXTRAN 1469 INTEGER IGTRAN(4,MXTRAN) 1470 INTEGER IGDOTS(MXVEC,MXTRAN) 1471 1472 LOGICAL LFNDC, LFNDA 1473 INTEGER IZETAV, ITAMPA, ITAMPB, ITAMPC 1474 INTEGER ITRAN, IVEC 1475 INTEGER ITAMP, I, IDX 1476 1477#include "priunit.h" 1478* statement functions: 1479 LOGICAL LGTST, LGEND 1480 INTEGER IL, IA, IB 1481 LGTST(ITRAN,IL,IA,IB) = IGTRAN(1,ITRAN).EQ.IL .AND. 1482 & (IGTRAN(2,ITRAN).EQ.IA .AND. IGTRAN(3,ITRAN).EQ.IB) 1483 LGEND(ITRAN) = ITRAN.GT.MXTRAN .OR. 1484 & (IGTRAN(1,ITRAN)+IGTRAN(2,ITRAN)+IGTRAN(3,ITRAN)).LE.0 1485 1486 1487 1488 IF (LOCDBG) THEN 1489 WRITE (LUPRI,*) 'CC_SETG212> on entry:' 1490 WRITE (LUPRI,*) 'IVEC :',IVEC 1491 WRITE (LUPRI,*) 'ITRAN:',ITRAN 1492 WRITE (LUPRI,*) 'ITAMPA,ITAMPB,ITAMPC:',ITAMPA,ITAMPB,ITAMPC 1493 IDX = 1 1494 DO WHILE (.NOT. LGEND(IDX)) 1495 WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') 'CC_SETG212>', 1496 & (IGTRAN(I,IDX),I=1,3),(IGDOTS(I,IDX),I=1,MXVEC) 1497 IDX = IDX + 1 1498 END DO 1499 CALL FLSHFO(LUPRI) 1500 END IF 1501 1502*---------------------------------------------------------------------* 1503* set up list of G matrix transformations 1504*---------------------------------------------------------------------* 1505 ITRAN = 1 1506 LFNDA = LGTST(ITRAN,IZETAV,ITAMPC,ITAMPB) 1507 LFNDC = LGTST(ITRAN,IZETAV,ITAMPA,ITAMPB) 1508 1509 DO WHILE ( .NOT. (LFNDA.OR.LFNDC.OR.LGEND(ITRAN)) ) 1510 ITRAN = ITRAN + 1 1511 LFNDA = LGTST(ITRAN,IZETAV,ITAMPC,ITAMPB) 1512 LFNDC = LGTST(ITRAN,IZETAV,ITAMPA,ITAMPB) 1513 END DO 1514 1515 IF (.NOT.(LFNDC.OR.LFNDA)) THEN 1516 IGTRAN(1,ITRAN) = IZETAV 1517 IGTRAN(2,ITRAN) = ITAMPA 1518 IGTRAN(3,ITRAN) = ITAMPB 1519 ITAMP = ITAMPC 1520 ELSE 1521 IF (LFNDC) ITAMP = ITAMPC 1522 IF (LFNDA) ITAMP = ITAMPA 1523 END IF 1524 1525 IVEC = 1 1526 DO WHILE (IGDOTS(IVEC,ITRAN).NE.ITAMP .AND. 1527 & IGDOTS(IVEC,ITRAN).NE.0 .AND. IVEC.LE.MXVEC) 1528 IVEC = IVEC + 1 1529 END DO 1530 1531 IGDOTS(IVEC,ITRAN) = ITAMP 1532 1533*---------------------------------------------------------------------* 1534 IF (LOCDBG .OR. IVEC.GT.MXVEC .OR. ITRAN.GT.MXTRAN) THEN 1535 WRITE (LUPRI,*) 'CC_SETG212> on exit:' 1536 WRITE (LUPRI,*) 'IVEC :',IVEC 1537 WRITE (LUPRI,*) 'ITRAN:',ITRAN 1538 WRITE (LUPRI,*) 'ITAMPA,ITAMPB,ITAMPC:',ITAMPA,ITAMPB,ITAMPC 1539 IDX = 1 1540 DO WHILE (.NOT. LGEND(IDX)) 1541 WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') 'CC_SETG212>', 1542 & (IGTRAN(I,IDX),I=1,3),(IGDOTS(I,IDX),I=1,MXVEC) 1543 IDX = IDX + 1 1544 END DO 1545 CALL FLSHFO(LUPRI) 1546 END IF 1547 1548 IF (IVEC.GT.MXVEC .OR. ITRAN.GT.MXTRAN) THEN 1549 CALL QUIT('Overflow error in CC_SETG212.') 1550 END IF 1551 1552 RETURN 1553 END 1554*---------------------------------------------------------------------* 1555* END OF SUBROUTINE CC_SETG212 * 1556*---------------------------------------------------------------------* 1557