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_4HYP */ 20*=====================================================================* 21 SUBROUTINE CC_4HYP(WORK,LWORK) 22*---------------------------------------------------------------------* 23* 24* Purpose: direct calculation of frequency-dependent fourth 25* hyperpolarizabilities for the Couple Cluster models 26* 27* CCS, CC2, CCSD 28* 29* Written by Christof Haettig, maj 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 "ccr1rsp.h" 40#include "cc5rinf.h" 41#include "cclists.h" 42 43* local parameters: 44 CHARACTER*(19) MSGDBG 45 PARAMETER (MSGDBG = '[debug] CC_4HYP> ') 46 LOGICAL LOCDBG 47 PARAMETER (LOCDBG = .FALSE. ) 48 49 CHARACTER*10 MODEL 50 LOGICAL LADD 51 INTEGER LWORK 52 INTEGER KHYPPOL, NBHYPPOL 53 INTEGER KEND0, LEND0, MXTRAN3, MXTRAN4, MXVEC2, MXVEC3 54 INTEGER K0HTRAN, K0HDOTS, K0HCONS, N0HTRAN, MX0HTRAN, MX0HDOTS 55 INTEGER K1HTRAN, K1HDOTS, K1HCONS, N1HTRAN, MX1HTRAN, MX1HDOTS 56 INTEGER K0GTRAN, K0GDOTS, K0GCONS, N0GTRAN, MX0GTRAN, MX0GDOTS 57 INTEGER K1GTRAN, K1GDOTS, K1GCONS, N1GTRAN, MX1GTRAN, MX1GDOTS 58 INTEGER K1FATRAN,K1FADOTS,K1FACONS,N1FATRAN,MX1FATRAN,MX1FADOTS 59 INTEGER KDTRAN, KDDOTS, KDCONS, NDTRAN, MXDTRAN, MXDDOTS 60 INTEGER KCTRAN, KCDOTS, KCCONS, NCTRAN, MXCTRAN, MXCDOTS 61 INTEGER KBTRAN, KBDOTS, KBCONS, NBTRAN, MXBTRAN, MXBDOTS 62 INTEGER KBATRAN, KBADOTS, KBACONS, NBATRAN, MXBATRAN, MXBADOTS 63 INTEGER K0FTRAN, K0FDOTS, K0FCONS, N0FTRAN, MX0FTRAN, MX0FDOTS 64 INTEGER KCHI3, KCHIDOTS,KCHICONS,NXTRAN, MXCHI3T, MXCHI3D 65 INTEGER 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 CLUST', 82 & 'ER PENTIC 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_4HYP called for unknown Coupled Cluster.') 95 END IF 96 97* print some debug/info output 98 IF (IPR5HYP .GT. 10) WRITE(LUPRI,*) 'CC_4HYP Workspace:',LWORK 99 100*---------------------------------------------------------------------* 101* allocate & initialize work space for hyperpolarizabilities 102*---------------------------------------------------------------------* 103 NBHYPPOL = N5ROPER * N5RFREQ 104 105 KHYPPOL = 1 106 KEND0 = KHYPPOL + 2*NBHYPPOL 107 108 109 MXVEC2 = (NLRTLBL+1) * NLRTLBL / 2 110 MXVEC3 = (NLRTLBL+2) * MXVEC2 / 3 111 MXTRAN3 = NLRTLBL * NLRTLBL * NLRTLBL 112 MXTRAN4 = NLRTLBL * NLRTLBL * NLRTLBL * NLRTLBL 113 MXVEC2 = MIN(MXVEC2,2*180*NBHYPPOL) ! 180 = # perm. for B^AB{O} 114 MXVEC3 = MIN(MXVEC3,2*180*NBHYPPOL) ! 180 = # perm. for B^AB{O} 115 MXTRAN3 = MIN(MXTRAN3,2*180*NBHYPPOL) ! 180 = # perm. for B^AB{O} 116 MXTRAN4 = MIN(MXTRAN4,2*180*NBHYPPOL) ! 180 = # perm. for B^AB{O} 117 118 MX0HTRAN = 5 * MXTRAN4 119 MX1HTRAN = 5 * MXTRAN4 120 MX0GTRAN = 4 * MXTRAN4 121 MX1GTRAN = 4 * MXTRAN4 122 MX1FATRAN = MXDIM_FATRAN * MXTRAN4 123 MXDTRAN = 5 * MXTRAN4 124 MXCTRAN = 4 * MXTRAN4 125 MXBTRAN = 3 * MXTRAN4 126 MXBATRAN = MXDIM_BATRAN * MXTRAN4 127 MX0FTRAN = 3 * MXTRAN3 128 MXCHI3T = 1 * MXTRAN3 129 130 MX0HDOTS = MXTRAN4 * MXVEC2 131 MX1HDOTS = MXTRAN4 * MXVEC2 132 MX0GDOTS = MXTRAN4 * MXVEC2 133 MX1GDOTS = MXTRAN4 * MXVEC2 134 MX1FADOTS = MXTRAN4 * MXVEC2 135 MXDDOTS = MXTRAN4 * MXVEC2 136 MXCDOTS = MXTRAN4 * MXVEC2 137 MXBDOTS = MXTRAN4 * MXVEC2 138 MXBADOTS = MXTRAN4 * MXVEC2 139 MX0FDOTS = MXTRAN3 * MXVEC3 140 MXCHI3D = MXTRAN3 * MXVEC3 141 142 143 K0HTRAN = KEND0 144 K0HDOTS = K0HTRAN + MX0HTRAN 145 K0HCONS = K0HDOTS + MX0HDOTS 146 KEND0 = K0HCONS + MX0HDOTS 147 148 K1HTRAN = KEND0 149 K1HDOTS = K1HTRAN + MX1HTRAN 150 K1HCONS = K1HDOTS + MX1HDOTS 151 KEND0 = K1HCONS + MX1HDOTS 152 153 K0GTRAN = KEND0 154 K0GDOTS = K0GTRAN + MX0GTRAN 155 K0GCONS = K0GDOTS + MX0GDOTS 156 KEND0 = K0GCONS + MX0GDOTS 157 158 K1GTRAN = KEND0 159 K1GDOTS = K1GTRAN + MX1GTRAN 160 K1GCONS = K1GDOTS + MX1GDOTS 161 KEND0 = K1GCONS + MX1GDOTS 162 163 K1FATRAN = KEND0 164 K1FADOTS = K1FATRAN + MX1FATRAN 165 K1FACONS = K1FADOTS + MX1FADOTS 166 KEND0 = K1FACONS + MX1FADOTS 167 168 KDTRAN = KEND0 169 KDDOTS = KDTRAN + MXDTRAN 170 KDCONS = KDDOTS + MXDDOTS 171 KEND0 = KDCONS + MXDDOTS 172 173 KCTRAN = KEND0 174 KCDOTS = KCTRAN + MXCTRAN 175 KCCONS = KCDOTS + MXCDOTS 176 KEND0 = KCCONS + MXCDOTS 177 178 KBTRAN = KEND0 179 KBDOTS = KBTRAN + MXBTRAN 180 KBCONS = KBDOTS + MXBDOTS 181 KEND0 = KBCONS + MXBDOTS 182 183 KBATRAN = KEND0 184 KBADOTS = KBATRAN + MXBATRAN 185 KBACONS = KBADOTS + MXBADOTS 186 KEND0 = KBACONS + MXBADOTS 187 188 K0FTRAN = KEND0 189 K0FDOTS = K0FTRAN + MX0FTRAN 190 K0FCONS = K0FDOTS + MX0FDOTS 191 KEND0 = K0FCONS + MX0FDOTS 192 193 KCHI3 = KEND0 194 KCHIDOTS = KCHI3 + MXCHI3T 195 KCHICONS = KCHIDOTS + MXCHI3D 196 KEND0 = KCHICONS + MXCHI3D 197 198 LEND0 = LWORK - KEND0 199 200 IF (LEND0.LT.0) THEN 201 WRITE (LUPRI,*) 'KEND0,LEND0:',KEND0,LEND0 202 CALL QUIT('Insufficient memory in CC_4HYP.') 203 END IF 204 205* initialize hyperpolarizabilities and the lists of contributions: 206 CALL DZERO(WORK(1),KEND0-1) 207 208 IF (LOCDBG) THEN 209 CALL DCOPY(MX0HDOTS,8.8888d88,0,WORK(K0HCONS),1) 210 CALL DCOPY(MX1HDOTS,8.8888d88,0,WORK(K1HCONS),1) 211 CALL DCOPY(MX0GDOTS,8.8888d88,0,WORK(K0GCONS),1) 212 CALL DCOPY(MX1GDOTS,8.8888d88,0,WORK(K1GCONS),1) 213 CALL DCOPY(MX1FADOTS,8.8888d88,0,WORK(K1FACONS),1) 214 CALL DCOPY(MXDDOTS,8.8888d88,0,WORK(KDCONS),1) 215 CALL DCOPY(MXCDOTS,8.8888d88,0,WORK(KCCONS),1) 216 CALL DCOPY(MXBDOTS,8.8888d88,0,WORK(KBCONS),1) 217 CALL DCOPY(MXBADOTS,8.8888d88,0,WORK(KBACONS),1) 218 CALL DCOPY(MX0FDOTS,8.8888d88,0,WORK(K0FCONS),1) 219 CALL DCOPY(MXCHI3D,8.8888d88,0,WORK(KCHICONS),1) 220 END IF 221 222*---------------------------------------------------------------------* 223* set up lists for H, G and F transformations and ETA vectors: 224*---------------------------------------------------------------------* 225 226 LADD = .FALSE. 227 228 CALL CC5R_SETUP(MXTRAN3, MXVEC3, MXTRAN4, MXVEC2, 229 & WORK(K0HTRAN), WORK(K0HDOTS), WORK(K0HCONS), N0HTRAN, 230 & WORK(K1HTRAN), WORK(K1HDOTS), WORK(K1HCONS), N1HTRAN, 231 & WORK(K0GTRAN), WORK(K0GDOTS), WORK(K0GCONS), N0GTRAN, 232 & WORK(K1GTRAN), WORK(K1GDOTS), WORK(K1GCONS), N1GTRAN, 233 & WORK(K1FATRAN),WORK(K1FADOTS),WORK(K1FACONS),N1FATRAN, 234 & WORK(KDTRAN), WORK(KDDOTS), WORK(KDCONS), NDTRAN, 235 & WORK(KCTRAN), WORK(KCDOTS), WORK(KCCONS), NCTRAN, 236 & WORK(KBTRAN), WORK(KBDOTS), WORK(KBCONS), NBTRAN, 237 & WORK(KBATRAN), WORK(KBADOTS), WORK(KBACONS), NBATRAN, 238 & WORK(K0FTRAN), WORK(K0FDOTS), WORK(K0FCONS), N0FTRAN, 239 & WORK(KCHI3), WORK(KCHIDOTS),WORK(KCHICONS),NXTRAN, 240 & WORK(KHYPPOL), NBHYPPOL, LADD ) 241 242*---------------------------------------------------------------------* 243* calculate H matrix contributions: 244*---------------------------------------------------------------------* 245 IOPT = 5 246 CALL CC_HMAT('L0 ','R2 ','R1 ','R1 ','R2 ',N0HTRAN, MXVEC2, 247 & WORK(K0HTRAN),WORK(K0HDOTS),WORK(K0HCONS), 248 & WORK(KEND0), LEND0, IOPT ) 249 250 IOPT = 5 251 CALL CC_HMAT('L1 ','R1 ','R1 ','R1 ','R2 ',N1HTRAN, MXVEC2, 252 & WORK(K1HTRAN),WORK(K1HDOTS),WORK(K1HCONS), 253 & WORK(KEND0), LEND0, IOPT ) 254 255*---------------------------------------------------------------------* 256* calculate G matrix contributions: 257*---------------------------------------------------------------------* 258 IOPT = 5 259 CALL CC_GMATRIX('L0 ','R2 ','R2 ','R2 ',N0GTRAN, MXVEC2, 260 & WORK(K0GTRAN),WORK(K0GDOTS),WORK(K0GCONS), 261 & WORK(KEND0), LEND0, IOPT ) 262 263 IOPT = 5 264 CALL CC_GMATRIX('L1 ','R2 ','R1 ','R2 ',N1GTRAN, MXVEC2, 265 & WORK(K1GTRAN),WORK(K1GDOTS),WORK(K1GCONS), 266 & WORK(KEND0), LEND0, IOPT ) 267 268*---------------------------------------------------------------------* 269* calculate F{O} matrix contributions: 270*---------------------------------------------------------------------* 271 CALL CCQR_FADRV('L1 ','o1 ','R2 ','R2 ',N1FATRAN, MXVEC2, 272 & WORK(K1FATRAN),WORK(K1FADOTS), 273 & WORK(K1FACONS), 274 & WORK(KEND0), LEND0, 'DOTP' ) 275 276*---------------------------------------------------------------------* 277* calculate D matrix contributions: 278*---------------------------------------------------------------------* 279 IOPT = 5 280 CALL CC_DMAT(WORK(KDTRAN),NDTRAN, 281 & 'R1 ','R1 ','R1 ','R1 ',IOPT,'L2 ', 282 & WORK(KDDOTS),WORK(KDCONS),MXVEC2,WORK(KEND0),LEND0) 283 284*---------------------------------------------------------------------* 285* calculate C matrix contributions: 286*---------------------------------------------------------------------* 287 IOPT = 5 288 CALL CC_CMAT(WORK(KCTRAN),NCTRAN,'R2 ','R1 ','R1 ',IOPT,'L2 ', 289 & WORK(KCDOTS),WORK(KCCONS),MXVEC2,WORK(KEND0),LEND0) 290 291*---------------------------------------------------------------------* 292* calculate B matrix contributions: 293*---------------------------------------------------------------------* 294 IOPT = 5 295 CALL CC_BMAT(WORK(KBTRAN),NBTRAN,'R2 ','R2 ',IOPT,'L2 ', 296 & WORK(KBDOTS),WORK(KBCONS),MXVEC2, 297 & .FALSE.,WORK(KEND0),LEND0) 298 299*---------------------------------------------------------------------* 300* calculate B{O} matrix contributions: 301*---------------------------------------------------------------------* 302 IOPT = 5 303 CALL CC_BAMAT(WORK(KBATRAN),NBATRAN,'o1 ','R2 ','R1 ',IOPT,'L2 ', 304 & WORK(KBADOTS),WORK(KBACONS),MXVEC2,WORK(KEND0),LEND0) 305 306*---------------------------------------------------------------------* 307* calculate F matrix contributions: 308*---------------------------------------------------------------------* 309 IOPT = 5 310 CALL CC_FMATRIX(WORK(K0FTRAN),N0FTRAN,'L0 ','R3 ',IOPT,'R3 ', 311 & WORK(K0FDOTS),WORK(K0FCONS),MXVEC3, 312 & WORK(KEND0), LEND0) 313 314*---------------------------------------------------------------------* 315* calculate chi3 vector contributions: 316*---------------------------------------------------------------------* 317 CALL CC_DOTDRV('X3 ','R3 ',NXTRAN,MXVEC3, 318 & WORK(KCHI3), WORK(KCHIDOTS), WORK(KCHICONS), 319 & WORK(KEND0), LEND0 ) 320 321*---------------------------------------------------------------------* 322* collect contributions and add them to hyperpolarizabilities 323*---------------------------------------------------------------------* 324 325 LADD = .TRUE. 326 327 CALL CC5R_SETUP(MXTRAN3, MXVEC3, MXTRAN4, MXVEC2, 328 & WORK(K0HTRAN), WORK(K0HDOTS), WORK(K0HCONS), N0HTRAN, 329 & WORK(K1HTRAN), WORK(K1HDOTS), WORK(K1HCONS), N1HTRAN, 330 & WORK(K0GTRAN), WORK(K0GDOTS), WORK(K0GCONS), N0GTRAN, 331 & WORK(K1GTRAN), WORK(K1GDOTS), WORK(K1GCONS), N1GTRAN, 332 & WORK(K1FATRAN),WORK(K1FADOTS),WORK(K1FACONS),N1FATRAN, 333 & WORK(KDTRAN), WORK(KDDOTS), WORK(KDCONS), NDTRAN, 334 & WORK(KCTRAN), WORK(KCDOTS), WORK(KCCONS), NCTRAN, 335 & WORK(KBTRAN), WORK(KBDOTS), WORK(KBCONS), NBTRAN, 336 & WORK(KBATRAN), WORK(KBADOTS), WORK(KBACONS), NBATRAN, 337 & WORK(K0FTRAN), WORK(K0FDOTS), WORK(K0FCONS), N0FTRAN, 338 & WORK(KCHI3), WORK(KCHIDOTS),WORK(KCHICONS),NXTRAN, 339 & WORK(KHYPPOL), NBHYPPOL, LADD ) 340 341*---------------------------------------------------------------------* 342* print output & return: 343*---------------------------------------------------------------------* 344 345 CALL CC5HYPPRT(WORK(KHYPPOL)) 346 347 RETURN 348 END 349 350*=====================================================================* 351* END OF SUBROUTINE CC_4HYP * 352*=====================================================================* 353c /* deck CC5HYPPRT */ 354*=====================================================================* 355 SUBROUTINE CC5HYPPRT(HYPERPOL) 356*---------------------------------------------------------------------* 357* 358* Purpose: print fourth hyperpolarizabilities 359* 360* Written by Christof Haettig in maj 1997. 361* 362*=====================================================================* 363#if defined (IMPLICIT_NONE) 364 IMPLICIT NONE 365#else 366# include "implicit.h" 367#endif 368#include "priunit.h" 369#include "ccorb.h" 370#include "ccsdinp.h" 371#include "cc5rinf.h" 372#include "ccroper.h" 373 374 375 CHARACTER*10 BLANKS 376 CHARACTER*80 STRING 377 INTEGER IFREQ, IOPER, ISYOLD, ISYTOT, IDX 378 379 380#if defined (SYS_CRAY) 381 REAL HYPERPOL(N5RFREQ,N5ROPER,2) 382 REAL HALF 383#else 384 DOUBLE PRECISION HYPERPOL(N5RFREQ,N5ROPER,2) 385 DOUBLE PRECISION HALF 386#endif 387 PARAMETER (HALF = 0.5D0) 388 389*---------------------------------------------------------------------* 390* print header for hyperpolarizability section 391*---------------------------------------------------------------------* 392 BLANKS = ' ' 393 STRING = ' RESULTS FOR THE FOURTH HYPERPOLARIZABILITIES ' 394 395 IF (CCS) THEN 396 CALL AROUND( BLANKS//'FINAL CCS'//STRING(1:45)//BLANKS ) 397 ELSE IF (CC2) THEN 398 CALL AROUND( BLANKS//'FINAL CC2'//STRING(1:45)//BLANKS ) 399 ELSE IF (CCSD) THEN 400 CALL AROUND( BLANKS//'FINAL CCSD'//STRING(1:45)//BLANKS ) 401 ELSE 402 CALL QUIT('CC5HYPPRT called for an unknown Coupled '// 403 & 'Cluster model.') 404 END IF 405 406 IF (IPR5HYP.GT.5) THEN 407 WRITE(LUPRI,'(/1X,6(1X,A," operator",5X),4X,A,8X,A,/,121("-"))') 408 & 'A','B','C','D','E','F','epsilon','(asy. Resp.)' 409 ELSE 410 WRITE(LUPRI,'(/1X,6(1X,A," operator",5X),4X,A,/,110("-"))') 411 & 'A','B','C','D','E','F','epsilon' 412 END IF 413 414 ISYOLD = 1 415 416 DO IOPER = 1, N5ROPER 417 ISYTOT = 1 418 419 DO IDX = 1, 6 420 ISYTOT = MULD2H(ISYTOT,ISYOPR(I5ROP(IOPER,IDX))) 421 END DO 422 423 IFREQ = 1 424 IF (ISYTOT.EQ.1) THEN 425 IF (IPR5HYP.GT.5) THEN 426 WRITE(LUPRI,'(/1X,6(A7,F8.4,1X),G16.8," (",G10.3,")")') 427 & (LBLOPR(I5ROP(IOPER,IDX))(1:7),FREQ5(IFREQ,IDX),IDX=1,6), 428 & HALF*(HYPERPOL(IFREQ,IOPER,1)+HYPERPOL(IFREQ,IOPER,2)), 429 & HALF*(HYPERPOL(IFREQ,IOPER,1)-HYPERPOL(IFREQ,IOPER,2)) 430 ELSE 431 WRITE(LUPRI,'(/1X,6(A7,F8.4,1X),G16.8)') 432 & (LBLOPR(I5ROP(IOPER,IDX))(1:7),FREQ5(IFREQ,IDX),IDX=1,6), 433 & HALF*(HYPERPOL(IFREQ,IOPER,1)+HYPERPOL(IFREQ,IOPER,2)) 434 END IF 435 ISYOLD = 1 436 ELSE 437 IF (IPR5HYP.GT.5) THEN 438 IF (ISYOLD.EQ.1) WRITE(LUPRI,*) 439 WRITE(LUPRI,'(1X,6(A7,A8,1X),6X,A,7X," (",5X,A,6X,")")') 440 & (LBLOPR(I5ROP(IOPER,IDX))(1:7),' -.- ',IDX=1,6), 441 & '---', 442 & '---' 443 ELSE IF (IPR5HYP.GT.0) THEN 444 IF (ISYOLD.EQ.1) WRITE(LUPRI,*) 445 WRITE(LUPRI,'(1X,6(A7,A8,1X),6X,A,7X)') 446 & (LBLOPR(I5ROP(IOPER,IDX))(1:7),' -.- ',IDX=1,6), 447 & '---' 448 END IF 449 ISYOLD = 0 450 END IF 451 452 DO IFREQ = 2, N5RFREQ 453 IF (ISYTOT.EQ.1) THEN 454 IF (IPR5HYP.GT.5) THEN 455 WRITE(LUPRI,'(1X,6(7X,F8.4,1X),G16.8," (",G10.3,")")') 456 & (FREQ5(IFREQ,IDX),IDX=1,6), 457 & HALF*(HYPERPOL(IFREQ,IOPER,1)+HYPERPOL(IFREQ,IOPER,2)), 458 & HALF*(HYPERPOL(IFREQ,IOPER,1)-HYPERPOL(IFREQ,IOPER,2)) 459 ELSE 460 WRITE(LUPRI,'(1X,6(7X,F8.4,1X),G16.8)') 461 & (FREQ5(IFREQ,IDX),IDX=1,6), 462 & HALF*(HYPERPOL(IFREQ,IOPER,1)+HYPERPOL(IFREQ,IOPER,2)) 463 END IF 464 END IF 465 END DO 466 467 END DO 468 469 RETURN 470 END 471*---------------------------------------------------------------------* 472* END OF SUBROUTINE CC5HYPPRT * 473*---------------------------------------------------------------------* 474c /* deck CC5R_SETUP */ 475*=====================================================================* 476 SUBROUTINE CC5R_SETUP(MXTRAN3, MXVEC3, MXTRAN4, MXVEC2, 477 & I0HTRAN, I0HDOTS, W0H, N0HTRAN, 478 & I1HTRAN, I1HDOTS, W1H, N1HTRAN, 479 & I0GTRAN, I0GDOTS, W0G, N0GTRAN, 480 & I1GTRAN, I1GDOTS, W1G, N1GTRAN, 481 & I1FATRAN, I1FADOTS, W1FA, N1FATRAN, 482 & IDTRAN, IDDOTS, WD, NDTRAN, 483 & ICTRAN, ICDOTS, WC, NCTRAN, 484 & IBTRAN, IBDOTS, WB, NBTRAN, 485 & IBATRAN, IBADOTS, WBA, NBATRAN, 486 & I0FTRAN, I0FDOTS, W0F, N0FTRAN, 487 & IXTRAN, IXDOTS, WCHI, NXTRAN, 488 & HYPPOL, MXHYPPOL, LADD ) 489*---------------------------------------------------------------------* 490* 491* Purpose: set up for CC5R section 492* - list of H^0 matrix transformations 493* - list of H^A matrix transformations 494* - list of G^0 matrix transformations 495* - list of G^A matrix transformations 496* - list of F^A{O} matrix transformations 497* - list of D matrix transformations 498* - list of C matrix transformations 499* - list of B matrix transformations 500* - list of B{O} matrix transformations 501* - list of F^0 matrix transformations 502* - list of chi3 dot products 503* 504* LADD = .FALSE. --> set lists 505* LADD = .TRUE. --> add contributions to hyperpolarizabilities 506* 507* Written by Christof Haettig, maj 1997. 508* 509*=====================================================================* 510#if defined (IMPLICIT_NONE) 511 IMPLICIT NONE 512#else 513# include "implicit.h" 514#endif 515#include "priunit.h" 516#include "ccorb.h" 517#include "cc5rinf.h" 518#include "cc5perm.h" 519#include "ccroper.h" 520#include "cclists.h" 521 522* local parameters: 523 CHARACTER*(20) MSGDBG 524 PARAMETER (MSGDBG = '[debug] CC5R_SETUP> ') 525 LOGICAL LOCDBG 526 PARAMETER (LOCDBG = .FALSE.) 527 528 INTEGER ORD1, ORD2, ORD3 529 PARAMETER ( ORD1 = 6 ) ! order (nb. of operators) 530 PARAMETER ( ORD2 = (ORD1-1)*ORD1/2 ) 531 PARAMETER ( ORD3 = (ORD1-2)*ORD2/3 ) 532 533 LOGICAL LADD 534 535 INTEGER IOP(ORD1), IR1(ORD1), IL1(ORD1), ISY(ORD1) 536 INTEGER IR2(ORD2), IL2(ORD2) 537 INTEGER IR3(ORD3), IX3(ORD3) 538 539 INTEGER MXVEC3, MXTRAN3, MXVEC2, MXTRAN4, MXHYPPOL 540 541 INTEGER I0HTRAN(5,MXTRAN4) 542 INTEGER I0HDOTS(MXVEC2,MXTRAN4) 543 INTEGER I1HTRAN(5,MXTRAN4) 544 INTEGER I1HDOTS(MXVEC2,MXTRAN4) 545 546 INTEGER I0GTRAN(4,MXTRAN4) 547 INTEGER I0GDOTS(MXVEC2,MXTRAN4) 548 INTEGER I1GTRAN(4,MXTRAN4) 549 INTEGER I1GDOTS(MXVEC2,MXTRAN4) 550 551 INTEGER I1FATRAN(MXDIM_FATRAN,MXTRAN4) 552 INTEGER I1FADOTS(MXVEC2,MXTRAN4) 553 554 INTEGER IDTRAN(5,MXTRAN4) 555 INTEGER IDDOTS(MXVEC2,MXTRAN4) 556 INTEGER ICTRAN(4,MXTRAN4) 557 INTEGER ICDOTS(MXVEC2,MXTRAN4) 558 INTEGER IBTRAN(3,MXTRAN4) 559 INTEGER IBDOTS(MXVEC2,MXTRAN4) 560 INTEGER IBATRAN(MXDIM_BATRAN,MXTRAN4) 561 INTEGER IBADOTS(MXVEC2,MXTRAN4) 562 563 INTEGER I0FTRAN(3,MXTRAN3) 564 INTEGER I0FDOTS(MXVEC3,MXTRAN3) 565 566 INTEGER IXTRAN(MXTRAN3) 567 INTEGER IXDOTS(MXVEC3,MXTRAN3) 568 569 INTEGER N0HTRAN, N0GTRAN, N0FTRAN 570 INTEGER N1HTRAN, N1GTRAN, N1FATRAN 571 INTEGER NDTRAN, NCTRAN, NBTRAN, NBATRAN, NXTRAN 572 INTEGER N5RRESF 573 INTEGER MXVH0, MXVH1, MXVG0, MXVG1, MXVF0, MXVFA1 574 INTEGER MXVD, MXVC, MXVB, MXVBA, MXCHI 575 576 INTEGER IVEC, ITRAN, I, IDX, IDXAB, IDXA, IDXB, IDXC, IDXABC 577 INTEGER ISYML, ISYM1, ISYM2, ISYM3, ISYTOT 578 INTEGER IFREQ, IOPER 579 INTEGER P, ISIGN 580 581#if defined (SYS_CRAY) 582 REAL SIGN, SUM, SSUM 583 REAL HYPPOL(MXHYPPOL) 584 REAL H0CON(45), H1CON(60) 585 REAL G0CON(15), G1CON(90) 586 REAL F0CON(10), F1ACON(90), CHICON(20) 587 REAL DCON(15), CCON(90), BCON(45), BACON(180) 588 REAL W0H(MXVEC2,MXTRAN4), W1H(MXVEC2,MXTRAN4) 589 REAL W0G(MXVEC2,MXTRAN4), W1G(MXVEC2,MXTRAN4) 590 REAL W1FA(MXVEC2,MXTRAN4), WD(MXVEC2,MXTRAN4) 591 REAL WC(MXVEC2,MXTRAN4), WB(MXVEC2,MXTRAN4) 592 REAL WBA(MXVEC2,MXTRAN4), W0F(MXVEC3,MXTRAN3) 593 REAL WCHI(MXVEC3,MXTRAN3) 594#else 595 DOUBLE PRECISION SIGN, SUM, SSUM 596 DOUBLE PRECISION HYPPOL(MXHYPPOL) 597 DOUBLE PRECISION H0CON(45), H1CON(60) 598 DOUBLE PRECISION G0CON(15), G1CON(90) 599 DOUBLE PRECISION F0CON(10), F1ACON(90), CHICON(20) 600 DOUBLE PRECISION DCON(15), CCON(90), BCON(45), BACON(180) 601 DOUBLE PRECISION W0H(MXVEC2,MXTRAN4), W1H(MXVEC2,MXTRAN4) 602 DOUBLE PRECISION W0G(MXVEC2,MXTRAN4), W1G(MXVEC2,MXTRAN4) 603 DOUBLE PRECISION W1FA(MXVEC2,MXTRAN4), WD(MXVEC2,MXTRAN4) 604 DOUBLE PRECISION WC(MXVEC2,MXTRAN4), WB(MXVEC2,MXTRAN4) 605 DOUBLE PRECISION WBA(MXVEC2,MXTRAN4), W0F(MXVEC3,MXTRAN3) 606 DOUBLE PRECISION WCHI(MXVEC3,MXTRAN3) 607#endif 608 609 610* external functions: 611 INTEGER IR3TAMP 612 INTEGER IR2TAMP 613 INTEGER IR1TAMP 614 INTEGER IL1ZETA 615 INTEGER IL2ZETA 616 INTEGER ICHI3 617 618 619*---------------------------------------------------------------------* 620* initializations: 621*---------------------------------------------------------------------* 622 IF (.NOT. LADD) THEN 623 N0HTRAN = 0 624 N1HTRAN = 0 625 N0GTRAN = 0 626 N1GTRAN = 0 627 N1FATRAN = 0 628 NDTRAN = 0 629 NCTRAN = 0 630 NBTRAN = 0 631 NBATRAN = 0 632 N0FTRAN = 0 633 NXTRAN = 0 634 END IF 635 636 MXVH0 = 0 637 MXVH1 = 0 638 MXVG0 = 0 639 MXVG1 = 0 640 MXVFA1 = 0 641 MXVD = 0 642 MXVC = 0 643 MXVB = 0 644 MXVBA = 0 645 MXVF0 = 0 646 MXCHI = 0 647 648 N5RRESF = 0 649 650 CALL DZERO(HYPPOL,MXHYPPOL) 651 652*---------------------------------------------------------------------* 653* start loop over all requested third hyperpolarizabilities 654*---------------------------------------------------------------------* 655 656 DO IOPER = 1, N5ROPER 657 ISYTOT = 1 658 DO IDX = 1, 6 659 IOP(IDX) = I5ROP(IOPER,IDX) 660 ISY(IDX) = ISYOPR(IOP(IDX)) 661 ISYTOT = MULD2H(ISYTOT,ISY(IDX)) 662 END DO 663 664 IF (LOCDBG) THEN 665 CALL DCOPY (45,999.99d0,0,H0CON,1) 666 CALL DCOPY (60,999.99d0,0,H1CON,1) 667 CALL DCOPY (15,999.99d0,0,G0CON,1) 668 CALL DCOPY (90,999.99d0,0,G1CON,1) 669 CALL DCOPY (10,999.99d0,0,F0CON,1) 670 CALL DCOPY (90,999.99d0,0,F1ACON,1) 671 CALL DCOPY (20,999.99d0,0,CHICON,1) 672 CALL DCOPY (60,999.99d0,0,DCON,1) 673 CALL DCOPY (90,999.99d0,0,CCON,1) 674 CALL DCOPY (45,999.99d0,0,BCON,1) 675 CALL DCOPY (180,999.99d0,0,BACON,1) 676 END IF 677 678 IF (ISYTOT .EQ. 1) THEN 679 680 DO IFREQ = 1, N5RFREQ 681 682 N5RRESF = N5RRESF + 1 683 684 DO ISIGN = 1, -1, -2 685 SIGN = DBLE(ISIGN) 686 687 DO IDX = 1, 6 688 IL1(IDX) = IL1ZETA(LBLOPR(IOP(IDX)),.FALSE., 689 & SIGN*FREQ5(IFREQ,IDX),ISYM1) 690 IR1(IDX) = IR1TAMP(LBLOPR(IOP(IDX)),.FALSE., 691 & SIGN*FREQ5(IFREQ,IDX),ISYM1) 692 END DO 693 694 IDXAB = 0 695 DO IDXB = 2, 6 696 DO IDXA = 1, IDXB-1 697 IDXAB = IDXAB + 1 698 IR2(IDXAB) = 699 & IR2TAMP(LBLOPR(IOP(IDXA)),.FALSE.,SIGN*FREQ5(IFREQ,IDXA),ISYM1, 700 & LBLOPR(IOP(IDXB)),.FALSE.,SIGN*FREQ5(IFREQ,IDXB),ISYM2) 701 702 IL2(IDXAB) = 703 & IL2ZETA(LBLOPR(IOP(IDXA)), SIGN*FREQ5(IFREQ,IDXA),ISYM1, 704 & LBLOPR(IOP(IDXB)), SIGN*FREQ5(IFREQ,IDXB),ISYM2) 705 END DO 706 END DO 707 708 IDXABC = 0 709 DO IDXC = 3, 6 710 DO IDXB = 2, IDXC-1 711 DO IDXA = 1, IDXB-1 712 IDXABC = IDXABC + 1 713 IR3(IDXABC) = 714 & IR3TAMP( LBLOPR(IOP(IDXA)),SIGN*FREQ5(IFREQ,IDXA),ISYM1, 715 & LBLOPR(IOP(IDXB)),SIGN*FREQ5(IFREQ,IDXB),ISYM2, 716 & LBLOPR(IOP(IDXC)),SIGN*FREQ5(IFREQ,IDXC),ISYM3) 717 718 IX3(IDXABC) = 719 & ICHI3( LBLOPR(IOP(IDXA)),SIGN*FREQ5(IFREQ,IDXA),ISYM1, 720 & LBLOPR(IOP(IDXB)),SIGN*FREQ5(IFREQ,IDXB),ISYM2, 721 & LBLOPR(IOP(IDXC)),SIGN*FREQ5(IFREQ,IDXC),ISYM3) 722 END DO 723 END DO 724 END DO 725 726*---------------------------------------------------------------------* 727* set up list of H^0 matrix transformations, 45 permutations 728*---------------------------------------------------------------------* 729 DO P = 1, 45 730 CALL CC_SETH2112(I0HTRAN,I0HDOTS,MXTRAN4,MXVEC2, 731 & 0,IR2(KP2(P)),IR1(K1(P)),IR1(K2(P)),IR2(KP3(P)),ITRAN,IVEC) 732 N0HTRAN = MAX(N0HTRAN,ITRAN) 733 MXVH0 = MAX(MXVH0,IVEC) 734 H0CON(P) = W0H(IVEC,ITRAN) 735 END DO 736 737*---------------------------------------------------------------------* 738* set up list of H^A matrix transformations, 60 permutations 739*---------------------------------------------------------------------* 740 DO P = 1, 15 741 CALL CC_SETH1112(I1HTRAN,I1HDOTS,MXTRAN4,MXVEC2, 742 & IL1(J1(P)),IR1(J2(P)),IR1(J3(P)),IR1(J4(P)),IR2(JP(P)), 743 & ITRAN,IVEC) 744 N1HTRAN = MAX(N1HTRAN,ITRAN) 745 MXVH1 = MAX(MXVH1,IVEC) 746 H1CON(P) = W1H(IVEC,ITRAN) 747 748 CALL CC_SETH1112(I1HTRAN,I1HDOTS,MXTRAN4,MXVEC2, 749 & IL1(J2(P)),IR1(J1(P)),IR1(J3(P)),IR1(J4(P)),IR2(JP(P)), 750 & ITRAN,IVEC) 751 N1HTRAN = MAX(N1HTRAN,ITRAN) 752 MXVH1 = MAX(MXVH1,IVEC) 753 H1CON(P+15) = W1H(IVEC,ITRAN) 754 755 CALL CC_SETH1112(I1HTRAN,I1HDOTS,MXTRAN4,MXVEC2, 756 & IL1(J3(P)),IR1(J1(P)),IR1(J2(P)),IR1(J4(P)),IR2(JP(P)), 757 & ITRAN,IVEC) 758 N1HTRAN = MAX(N1HTRAN,ITRAN) 759 MXVH1 = MAX(MXVH1,IVEC) 760 H1CON(P+30) = W1H(IVEC,ITRAN) 761 762 CALL CC_SETH1112(I1HTRAN,I1HDOTS,MXTRAN4,MXVEC2, 763 & IL1(J4(P)),IR1(J1(P)),IR1(J2(P)),IR1(J3(P)),IR2(JP(P)), 764 & ITRAN,IVEC) 765 N1HTRAN = MAX(N1HTRAN,ITRAN) 766 MXVH1 = MAX(MXVH1,IVEC) 767 H1CON(P+45) = W1H(IVEC,ITRAN) 768 END DO 769 770*---------------------------------------------------------------------* 771* set up list of G^0 matrix transformations, 15 permutations 772*---------------------------------------------------------------------* 773 DO P = 1, 15 774 CALL CCQR_SETG(I0GTRAN,I0GDOTS,MXTRAN4,MXVEC2, 775 & 0,IR2(IP1(P)),IR2(IP2(P)),IR2(IP3(P)),ITRAN,IVEC) 776 N0GTRAN = MAX(N0GTRAN,ITRAN) 777 MXVG0 = MAX(MXVG0,IVEC) 778 G0CON(P) = W0G(IVEC,ITRAN) 779 END DO 780 781*---------------------------------------------------------------------* 782* set up list of G^A matrix transformations, 90 permutations 783*---------------------------------------------------------------------* 784 DO P = 1, 45 785 CALL CC_SETG212(I1GTRAN,I1GDOTS,MXTRAN4,MXVEC2, 786 & IL1(K1(P)),IR2(KP2(P)),IR1(K2(P)),IR2(KP3(P)),ITRAN,IVEC) 787 N1GTRAN = MAX(N1GTRAN,ITRAN) 788 MXVG1 = MAX(MXVG1,IVEC) 789 G1CON(P) = W1G(IVEC,ITRAN) 790 791 CALL CC_SETG212(I1GTRAN,I1GDOTS,MXTRAN4,MXVEC2, 792 & IL1(K2(P)),IR2(KP2(P)),IR1(K1(P)),IR2(KP3(P)),ITRAN,IVEC) 793 N1GTRAN = MAX(N1GTRAN,ITRAN) 794 MXVG1 = MAX(MXVG1,IVEC) 795 G1CON(P+45) = W1G(IVEC,ITRAN) 796 END DO 797 798*---------------------------------------------------------------------* 799* set up list of F^A{O} matrix transformations, 90 permutations 800*---------------------------------------------------------------------* 801 DO P = 1, 45 802 CALL CCQR_SETFA(I1FATRAN,I1FADOTS,MXTRAN4,MXVEC2, 803 & IL1(K1(P)),IOP(K2(P)),IR2(KP2(P)),IR2(KP3(P)),ITRAN,IVEC) 804 N1FATRAN = MAX(N1FATRAN,ITRAN) 805 MXVFA1 = MAX(MXVFA1,IVEC) 806 F1ACON(P) = W1FA(IVEC,ITRAN) 807 808 CALL CCQR_SETFA(I1FATRAN,I1FADOTS,MXTRAN4,MXVEC2, 809 & IL1(K2(P)),IOP(K1(P)),IR2(KP2(P)),IR2(KP3(P)),ITRAN,IVEC) 810 N1FATRAN = MAX(N1FATRAN,ITRAN) 811 MXVFA1 = MAX(MXVFA1,IVEC) 812 F1ACON(P+45) = W1FA(IVEC,ITRAN) 813 END DO 814 815*---------------------------------------------------------------------* 816* set up list of D matrix transformations, 15 permutations 817*---------------------------------------------------------------------* 818 DO P = 1, 15 819 CALL CC_SETD1111(IDTRAN,IDDOTS,MXTRAN4,MXVEC2, 820 & IL2(JP(P)),IR1(J1(P)),IR1(J2(P)),IR1(J3(P)),IR1(J4(P)), 821 & ITRAN,IVEC) 822 NDTRAN = MAX(NDTRAN,ITRAN) 823 MXVD = MAX(MXVD,IVEC) 824 DCON(P) = WD(IVEC,ITRAN) 825 END DO 826 827*---------------------------------------------------------------------* 828* set up list of C matrix transformations, 90 permutations 829*---------------------------------------------------------------------* 830 DO P = 1, 45 831 CALL CC_SETC211(ICTRAN,ICDOTS,MXTRAN4,MXVEC2, 832 & IL2(KP2(P)),IR2(KP3(P)),IR1(K1(P)),IR1(K2(P)),ITRAN,IVEC) 833 NCTRAN = MAX(NCTRAN,ITRAN) 834 MXVC = MAX(MXVC,IVEC) 835 CCON(P) = WC(IVEC,ITRAN) 836 837 CALL CC_SETC211(ICTRAN,ICDOTS,MXTRAN4,MXVEC2, 838 & IL2(KP3(P)),IR2(KP2(P)),IR1(K1(P)),IR1(K2(P)),ITRAN,IVEC) 839 NCTRAN = MAX(NCTRAN,ITRAN) 840 MXVC = MAX(MXVC,IVEC) 841 CCON(P+45) = WC(IVEC,ITRAN) 842 END DO 843 844*---------------------------------------------------------------------* 845* set up list of B matrix transformations, 45 permutations 846*---------------------------------------------------------------------* 847 DO P = 1, 45 848 CALL CC_SETB11(IBTRAN,IBDOTS,MXTRAN4,MXVEC2, 849 & IL2(KP1(P)),IR2(KP2(P)),IR2(KP3(P)),ITRAN,IVEC) 850 NBTRAN = MAX(NBTRAN,ITRAN) 851 MXVB = MAX(MXVB,IVEC) 852 BCON(P) = WB(IVEC,ITRAN) 853 END DO 854 855*---------------------------------------------------------------------* 856* set up list of B{O} matrix transformations, 180 permutations 857*---------------------------------------------------------------------* 858 DO P = 1, 45 859 CALL CC_SETBA12(IBATRAN,IBADOTS,MXTRAN4,MXVEC2, 860 & IL2(KP2(P)),IOP(K1(P)),IR2(KP3(P)),IR1(K2(P)),ITRAN,IVEC) 861 NBATRAN = MAX(NBATRAN,ITRAN) 862 MXVBA = MAX(MXVBA,IVEC) 863 BACON(P) = WBA(IVEC,ITRAN) 864 865 CALL CC_SETBA12(IBATRAN,IBADOTS,MXTRAN4,MXVEC2, 866 & IL2(KP3(P)),IOP(K1(P)),IR2(KP2(P)),IR1(K2(P)),ITRAN,IVEC) 867 NBATRAN = MAX(NBATRAN,ITRAN) 868 MXVBA = MAX(MXVBA,IVEC) 869 BACON(P+45) = WBA(IVEC,ITRAN) 870 871 CALL CC_SETBA12(IBATRAN,IBADOTS,MXTRAN4,MXVEC2, 872 & IL2(KP2(P)),IOP(K2(P)),IR2(KP3(P)),IR1(K1(P)),ITRAN,IVEC) 873 NBATRAN = MAX(NBATRAN,ITRAN) 874 MXVBA = MAX(MXVBA,IVEC) 875 BACON(P+90) = WBA(IVEC,ITRAN) 876 877 CALL CC_SETBA12(IBATRAN,IBADOTS,MXTRAN4,MXVEC2, 878 & IL2(KP3(P)),IOP(K2(P)),IR2(KP2(P)),IR1(K1(P)),ITRAN,IVEC) 879 NBATRAN = MAX(NBATRAN,ITRAN) 880 MXVBA = MAX(MXVBA,IVEC) 881 BACON(P+135) = WBA(IVEC,ITRAN) 882 END DO 883 884*---------------------------------------------------------------------* 885* set up list of F^0 matrix transformations, 10 permutations 886*---------------------------------------------------------------------* 887 DO P = 1, 10 888 CALL CCQR_SETF(I0FTRAN,I0FDOTS,MXTRAN3,MXVEC3, 889 & 0,IR3(IT1(P)),IR3(IT2(P)),ITRAN,IVEC) 890 N0FTRAN = MAX(N0FTRAN,ITRAN) 891 MXVF0 = MAX(MXVF0,IVEC) 892 F0CON(P) = W0F(IVEC,ITRAN) 893 END DO 894 895*---------------------------------------------------------------------* 896* set up list of chi3 vector dot products, 20 permutations 897*---------------------------------------------------------------------* 898 DO P = 1, 10 899 CALL CC_SETDOT(IXTRAN,IXDOTS,MXTRAN3,MXVEC3, 900 & IX3(IT1(P)),IR3(IT2(P)),ITRAN,IVEC) 901 NXTRAN = MAX(NXTRAN,ITRAN) 902 MXCHI = MAX(MXCHI,IVEC) 903 CHICON(P) = WCHI(IVEC,ITRAN) 904 905 CALL CC_SETDOT(IXTRAN,IXDOTS,MXTRAN3,MXVEC3, 906 & IX3(IT2(P)),IR3(IT1(P)),ITRAN,IVEC) 907 NXTRAN = MAX(NXTRAN,ITRAN) 908 MXCHI = MAX(MXCHI,IVEC) 909 CHICON(P+10) = WCHI(IVEC,ITRAN) 910 END DO 911 912*---------------------------------------------------------------------* 913* add all 660 contributions to the hyperpolarizabilities: 914*---------------------------------------------------------------------* 915 IF (LADD) THEN 916 IDX = N5RFREQ*(IOPER-1) + IFREQ 917 IF (ISIGN.EQ.-1) IDX = IDX + N5RFREQ*N5ROPER 918 919 DO P = 1, 10 920 HYPPOL(IDX) = HYPPOL(IDX) + F0CON(P) + CHICON(P)+CHICON(P+10) 921 END DO 922 923 DO P = 1, 15 924 HYPPOL(IDX) = HYPPOL(IDX) + G0CON(P) + DCON(P) + 925 & H1CON(P) + H1CON(P+15) + H1CON(P+30) + H1CON(P+45) 926 END DO 927 928 DO P = 1, 45 929 HYPPOL(IDX) = HYPPOL(IDX) + 930 & H0CON(P) + G1CON(P) + G1CON(P+45) + 931 & F1ACON(P) + F1ACON(P+45) + CCON(P) + CCON(P+45) + BCON(P) + 932 & BACON(P) + BACON(P+45) + BACON(P+90) + BACON(P+135) 933 END DO 934 935 IF (LOCDBG) THEN 936 WRITE (LUPRI,*) 'IOP:',IOP 937 WRITE (LUPRI,*) 'IDX',IDX 938 WRITE (LUPRI,*) 'HYPPOL:',HYPPOL(IDX) 939 940 SSUM= 0.0d0 941 942 SUM = 0.0d0 943 DO I = 1, 10 944 SUM = SUM + F0CON(I) 945 END DO 946 SSUM = SSUM + SUM 947 WRITE (LUPRI,*) 'F0CON:',F0CON 948 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 949 950 SUM = 0.0d0 951 DO I = 1, 20 952 SUM = SUM + CHICON(I) 953 END DO 954 SSUM = SSUM + SUM 955 WRITE (LUPRI,*) 'CHICON:',CHICON 956 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 957 958 SUM = 0.0d0 959 DO I = 1, 15 960 SUM = SUM + G0CON(I) 961 END DO 962 SSUM = SSUM + SUM 963 WRITE (LUPRI,*) 'G0CON:',G0CON 964 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 965 966 SUM = 0.0d0 967 DO I = 1, 15 968 SUM = SUM + DCON(I) 969 END DO 970 SSUM = SSUM + SUM 971 WRITE (LUPRI,*) 'DCON:',DCON 972 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 973 974 SUM = 0.0d0 975 DO I = 1, 60 976 SUM = SUM + H1CON(I) 977 END DO 978 SSUM = SSUM + SUM 979 WRITE (LUPRI,*) 'H1CON:',H1CON 980 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 981 982 SUM = 0.0d0 983 DO I = 1, 45 984 SUM = SUM + H0CON(I) 985 END DO 986 SSUM = SSUM + SUM 987 WRITE (LUPRI,*) 'H0CON:',H0CON 988 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 989 990 SUM = 0.0d0 991 DO I = 1, 90 992 SUM = SUM + G1CON(I) 993 END DO 994 SSUM = SSUM + SUM 995 WRITE (LUPRI,*) 'G1CON:',G1CON 996 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 997 998 SUM = 0.0d0 999 DO I = 1, 90 1000 SUM = SUM + F1ACON(I) 1001 END DO 1002 SSUM = SSUM + SUM 1003 WRITE (LUPRI,*) 'F1ACON:',F1ACON 1004 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1005 1006 SUM = 0.0d0 1007 DO I = 1, 90 1008 SUM = SUM + CCON(I) 1009 END DO 1010 SSUM = SSUM + SUM 1011 WRITE (LUPRI,*) 'CCON:',CCON 1012 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1013 1014 SUM = 0.0d0 1015 DO I = 1, 45 1016 SUM = SUM + BCON(I) 1017 END DO 1018 SSUM = SSUM + SUM 1019 WRITE (LUPRI,*) 'BCON:',BCON 1020 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1021 1022 SUM = 0.0d0 1023 DO I = 1, 180 1024 SUM = SUM + BACON(I) 1025 END DO 1026 SSUM = SSUM + SUM 1027 WRITE (LUPRI,*) 'BACON:',BACON 1028 WRITE (LUPRI,*) 'SUM:', SUM, SSUM 1029 1030 WRITE (LUPRI,*) 'HYPPOL:',HYPPOL(IDX) 1031 END IF 1032 1033 END IF 1034 1035*---------------------------------------------------------------------* 1036* end loop over all requested hyperpolarizabilities 1037*---------------------------------------------------------------------* 1038 END DO 1039 END DO 1040 END IF 1041 END DO 1042 1043*---------------------------------------------------------------------* 1044* print statistics and lists: 1045*---------------------------------------------------------------------* 1046 IF (.NOT. LADD) THEN 1047 1048* general statistics: 1049 WRITE(LUPRI,'(/,/3X,A,I4,A)') 'For the requested',N5RRESF, 1050 & ' pentic response functions ' 1051 WRITE(LUPRI,'((8X,A,I3,A))') 1052 & ' - ',N0HTRAN, ' H^0 matrix transformations ', 1053 & ' - ',N1HTRAN, ' H^A matrix transformations ', 1054 & ' - ',N0GTRAN, ' G^0 matrix transformations ', 1055 & ' - ',N1GTRAN, ' G^A matrix transformations ', 1056 & ' - ',N1FATRAN, ' F^A{O} matrix transformations ', 1057 & ' - ',NDTRAN, ' D matrix transformations ', 1058 & ' - ',NCTRAN, ' C matrix transformations ', 1059 & ' - ',NBTRAN, ' B matrix transformations ', 1060 & ' - ',NBATRAN, ' B{O} matrix transformations ', 1061 & ' - ',N0FTRAN, ' F^0 matrix transformations ', 1062 & ' - ',NXTRAN, ' dot products with chi3 vectors' 1063 WRITE(LUPRI,'(3X,A,/,/)') 'will be performed.' 1064 1065 1066* H^0 matrix transformations: 1067 IF (LOCDBG) THEN 1068 WRITE (LUPRI,*) 'List of H^0 matrix transformations:' 1069 DO ITRAN = 1, N0HTRAN 1070 WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG, 1071 & (I0HTRAN(I,ITRAN),I=1,4),(I0HDOTS(I,ITRAN),I=1,MXVH0) 1072 END DO 1073 WRITE (LUPRI,*) 1074 END IF 1075 1076* H^A matrix transformations: 1077 IF (LOCDBG) THEN 1078 WRITE (LUPRI,*) 'List of H^A matrix transformations:' 1079 DO ITRAN = 1, N1HTRAN 1080 WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG, 1081 & (I1HTRAN(I,ITRAN),I=1,4),(I1HDOTS(I,ITRAN),I=1,MXVH1) 1082 END DO 1083 WRITE (LUPRI,*) 1084 END IF 1085 1086* G^0 matrix transformations: 1087 IF (LOCDBG) THEN 1088 WRITE (LUPRI,*) 'List of G^0 matrix transformations:' 1089 DO ITRAN = 1, N0GTRAN 1090 WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG, 1091 & (I0GTRAN(I,ITRAN),I=1,3),(I0GDOTS(I,ITRAN),I=1,MXVG0) 1092 END DO 1093 WRITE (LUPRI,*) 1094 END IF 1095 1096* G^A matrix transformations: 1097 IF (LOCDBG) THEN 1098 WRITE (LUPRI,*) 'List of G^A matrix transformations:' 1099 DO ITRAN = 1, N1GTRAN 1100 WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG, 1101 & (I1GTRAN(I,ITRAN),I=1,3),(I1GDOTS(I,ITRAN),I=1,MXVG1) 1102 END DO 1103 WRITE (LUPRI,*) 1104 END IF 1105 1106* F^0 matrix transformations: 1107 IF (LOCDBG) THEN 1108 WRITE (LUPRI,*) 'List of F^0 matrix transformations:' 1109 DO ITRAN = 1, N0FTRAN 1110 WRITE(LUPRI,'(A,2I5,5X,(12I5,20X))') MSGDBG, 1111 & (I0FTRAN(I,ITRAN),I=1,2),(I0FDOTS(I,ITRAN),I=1,MXVF0) 1112 END DO 1113 WRITE (LUPRI,*) 1114 END IF 1115 1116* F^A{O} matrix transformations: 1117 IF (LOCDBG) THEN 1118 WRITE (LUPRI,*) 'List of F^A{O} matrix transformations:' 1119 DO ITRAN = 1, N1FATRAN 1120 WRITE(LUPRI,'(A,5I5,5X,(12I5,20X))') MSGDBG, 1121 & (I1FATRAN(I,ITRAN),I=1,5),(I1FADOTS(I,ITRAN),I=1,MXVFA1) 1122 END DO 1123 WRITE (LUPRI,*) 1124 END IF 1125 1126* D matrix transformations: 1127 IF (LOCDBG) THEN 1128 WRITE (LUPRI,*) 'List of D matrix transformations:' 1129 DO ITRAN = 1, NDTRAN 1130 WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG, 1131 & (IDTRAN(I,ITRAN),I=1,4),(IDDOTS(I,ITRAN),I=1,MXVD) 1132 END DO 1133 WRITE (LUPRI,*) 1134 END IF 1135 1136* C matrix transformations: 1137 IF (LOCDBG) THEN 1138 WRITE (LUPRI,*) 'List of C matrix transformations:' 1139 DO ITRAN = 1, NCTRAN 1140 WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG, 1141 & (ICTRAN(I,ITRAN),I=1,3),(ICDOTS(I,ITRAN),I=1,MXVC) 1142 END DO 1143 WRITE (LUPRI,*) 1144 END IF 1145 1146* B matrix transformations: 1147 IF (LOCDBG) THEN 1148 WRITE (LUPRI,*) 'List of B matrix transformations:' 1149 DO ITRAN = 1, NBTRAN 1150 WRITE(LUPRI,'(A,2I5,5X,(12I5,20X))') MSGDBG, 1151 & (IBTRAN(I,ITRAN),I=1,2),(IBDOTS(I,ITRAN),I=1,MXVB) 1152 END DO 1153 WRITE (LUPRI,*) 1154 END IF 1155 1156* B{O} matrix transformations: 1157 IF (LOCDBG) THEN 1158 WRITE (LUPRI,*) 'List of B{A} matrix transformations:' 1159 DO ITRAN = 1, NBATRAN 1160 WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG, 1161 & (IBATRAN(I,ITRAN),I=1,3),(IBADOTS(I,ITRAN),I=1,MXVBA) 1162 END DO 1163 WRITE (LUPRI,*) 1164 END IF 1165 1166* chi vector dot products: 1167 IF (LOCDBG) THEN 1168 WRITE (LUPRI,*) 'List of chi vector dot products:' 1169 DO ITRAN = 1, NXTRAN 1170 WRITE(LUPRI,'(A,I5,5X,(12I5,20X))') MSGDBG, 1171 & IXTRAN(ITRAN),(IXDOTS(I,ITRAN),I=1,MXCHI) 1172 END DO 1173 WRITE (LUPRI,*) 1174 END IF 1175 1176 END IF ! (.NOT. LADD) 1177*---------------------------------------------------------------------* 1178* return 1179*---------------------------------------------------------------------* 1180 1181 RETURN 1182 END 1183 1184*---------------------------------------------------------------------* 1185* END OF SUBROUTINE CC5R_SETUP * 1186*---------------------------------------------------------------------* 1187c /* deck CC_SETBA12 */ 1188*=====================================================================* 1189 SUBROUTINE CC_SETBA12(IBATRAN,IBADOTS,MXTRAN,MXVEC, 1190 & IZETAV,IOPER,ITAMPA,ITAMPB,ITRAN,IVEC) 1191*---------------------------------------------------------------------* 1192* 1193* Purpose: maintains a list of dot products of B{O} matrix 1194* transformations with two right amplitude vectors: 1195* L * (B{O} * T^A * T^B) 1196* assumes that T^A and T^B belong to different lists 1197* 1198* IBATRAN - list of F matrix transformations 1199* IBADOTS - list of vectors it should be dottet on 1200* 1201* MXTRAN - maximum list dimension 1202* MXVEC - maximum second dimension for IBADOTS 1203* 1204* IZETAV - index of lagrangian multiplier vector 1205* IOPER - index of property operator 1206* ITAMPA - index of amplitude vector A 1207* ITAMPB - index of amplitude vector B 1208* 1209* ITRAN - index in IBATRAN list 1210* IVEC - second index in IBADOTS list 1211* 1212* Written by Christof Haettig, june 1997. 1213* 1214*=====================================================================* 1215 IMPLICIT NONE 1216#include "priunit.h" 1217#include "cclists.h" 1218 1219 INTEGER MXVEC, MXTRAN 1220 INTEGER IBATRAN(MXDIM_BATRAN,MXTRAN) 1221 INTEGER IBADOTS(MXVEC,MXTRAN) 1222 1223 LOGICAL LFNDB 1224 INTEGER IZETAV, IOPER, ITAMPA, ITAMPB 1225 INTEGER ITRAN, IVEC 1226 INTEGER ITAMP, I, IDX 1227 1228* statement functions: 1229 LOGICAL LBATST, LBAEND 1230 INTEGER IB, IA, IO 1231 LBATST(ITRAN,IO,IA,IB) = IBATRAN(1,ITRAN).EQ.IO 1232 & .AND. IBATRAN(2,ITRAN).EQ.IA .AND. IBATRAN(3,ITRAN).EQ.IB 1233 LBAEND(ITRAN) = ITRAN.GT.MXTRAN .OR. 1234 & (IBATRAN(1,ITRAN)+IBATRAN(2,ITRAN)+IBATRAN(3,ITRAN)).LE.0 1235 1236 1237*---------------------------------------------------------------------* 1238* set up list of B{A} matrix transformations 1239*---------------------------------------------------------------------* 1240 ITRAN = 1 1241 LFNDB = LBATST(ITRAN,IOPER,ITAMPA,ITAMPB) 1242 1243 DO WHILE ( .NOT. (LFNDB.OR.LBAEND(ITRAN))) 1244 ITRAN = ITRAN + 1 1245 LFNDB = LBATST(ITRAN,IOPER,ITAMPA,ITAMPB) 1246 END DO 1247 1248 IF (.NOT.LFNDB) THEN 1249 IBATRAN(1,ITRAN) = IOPER 1250 IBATRAN(2,ITRAN) = ITAMPA 1251 IBATRAN(3,ITRAN) = ITAMPB 1252 IBATRAN(4,ITRAN) = 0 1253 IBATRAN(5,ITRAN) = 0 1254 IBATRAN(6,ITRAN) = 0 1255 IBATRAN(7,ITRAN) = 0 1256 IBATRAN(8,ITRAN) = 0 1257 ITAMP = IZETAV 1258 ELSE 1259 IF (LFNDB) ITAMP = IZETAV 1260 END IF 1261 1262 IVEC = 1 1263 DO WHILE (IBADOTS(IVEC,ITRAN).NE.ITAMP .AND. 1264 & IBADOTS(IVEC,ITRAN).NE.0 .AND. IVEC.LE.MXVEC) 1265 IVEC = IVEC + 1 1266 END DO 1267 1268 IBADOTS(IVEC,ITRAN) = ITAMP 1269 1270C WRITE (LUPRI,*) 'CC_SETBA12>',IZETAV,IOPER,ITAMPA,ITAMPB,ITRAN,IVEC 1271*---------------------------------------------------------------------* 1272 IF (IVEC.GT.MXVEC .OR. ITRAN.GT.MXTRAN) THEN 1273 WRITE (LUPRI,*) 'IZETAV,IOPER,ITAMPA,ITAMPB,ITRAN,IVEC:', 1274 & IZETAV,IOPER,ITAMPA,ITAMPB,ITRAN,IVEC 1275 IDX = 1 1276 DO WHILE ( .NOT. LBAEND(IDX) ) 1277 WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') 'CC_SETBA12>', 1278 & (IBATRAN(I,IDX),I=1,4),(IBADOTS(I,IDX),I=1,MXVEC) 1279 IDX = IDX + 1 1280 END DO 1281 CALL FLSHFO(LUPRI) 1282 CALL QUIT('Overflow error in CC_SETBA12.') 1283 END IF 1284 1285 RETURN 1286 END 1287 1288*---------------------------------------------------------------------* 1289* END OF SUBROUTINE CC_SETBA12 * 1290*---------------------------------------------------------------------* 1291c /* deck CC_SETB11 */ 1292*=====================================================================* 1293 SUBROUTINE CC_SETB11(IBTRAN,IBDOTS,MXTRAN,MXVEC, 1294 & IZETAV,ITAMPA,ITAMPB,ITRAN,IVEC) 1295*---------------------------------------------------------------------* 1296* 1297* Purpose: set up list of B matrix transformations 1298* assumes that T^A and T^B members of the same list 1299* 1300* IBTRAN - list of B matrix transformations 1301* IBDOTS - list of vectors it should be dottet on 1302* 1303* MXTRAN - maximum list dimension 1304* MXVEC - maximum second dimension for IBDOTS 1305* 1306* IZETAV - index of lagrangian multiplier vector 1307* ITAMPA - index of amplitude vector A 1308* ITAMPB - index of amplitude vector B 1309* 1310* ITRAN - index in IBTRAN list 1311* IVEC - second index in IBDOTS list 1312* 1313* Written by Christof Haettig, june 1997. 1314* 1315*=====================================================================* 1316 IMPLICIT NONE 1317 1318#include "priunit.h" 1319 INTEGER MXVEC, MXTRAN 1320 INTEGER IBTRAN(3,MXTRAN) 1321 INTEGER IBDOTS(MXVEC,MXTRAN) 1322 1323 LOGICAL LFND 1324 INTEGER IZETAV, ITAMPA, ITAMPB 1325 INTEGER ITRAN, IVEC 1326 INTEGER ITAMP, I, IDX 1327 1328* statement functions: 1329 LOGICAL LBTST, LBEND 1330 INTEGER IB, IA 1331 LBTST(ITRAN,IA,IB) = 1332 & ( IBTRAN(1,ITRAN).EQ.IA .AND. IBTRAN(2,ITRAN).EQ.IB ) .OR. 1333 & ( IBTRAN(1,ITRAN).EQ.IB .AND. IBTRAN(2,ITRAN).EQ.IA ) 1334 LBEND(ITRAN) = ITRAN.GT.MXTRAN .OR. 1335 & (IBTRAN(1,ITRAN)+IBTRAN(2,ITRAN)).LE.0 1336 1337*---------------------------------------------------------------------* 1338* set up list of B matrix transformations 1339*---------------------------------------------------------------------* 1340 ITRAN = 1 1341 LFND = LBTST(ITRAN,ITAMPA,ITAMPB) 1342 1343 DO WHILE ( .NOT. (LFND.OR.LBEND(ITRAN)) ) 1344 ITRAN = ITRAN + 1 1345 LFND = LBTST(ITRAN,ITAMPA,ITAMPB) 1346 END DO 1347 1348 IF (.NOT.LFND) THEN 1349 IBTRAN(1,ITRAN) = ITAMPA 1350 IBTRAN(2,ITRAN) = ITAMPB 1351 IBTRAN(3,ITRAN) = 0 1352 ITAMP = IZETAV 1353 ELSE 1354 IF (LFND) ITAMP = IZETAV 1355 END IF 1356 1357 IVEC = 1 1358 DO WHILE (IBDOTS(IVEC,ITRAN).NE.ITAMP .AND. 1359 & IBDOTS(IVEC,ITRAN).NE.0 .AND. IVEC.LE.MXVEC) 1360 IVEC = IVEC + 1 1361 END DO 1362 1363 IBDOTS(IVEC,ITRAN) = ITAMP 1364 1365*---------------------------------------------------------------------* 1366 IF (IVEC.GT.MXVEC .OR. ITRAN.GT.MXTRAN) THEN 1367 WRITE (LUPRI,*) 'IVEC :',IVEC 1368 WRITE (LUPRI,*) 'ITRAN:',ITRAN 1369 WRITE (LUPRI,*) 'ITAMPA,ITAMPB:',ITAMPA,ITAMPB 1370 IDX = 1 1371 DO WHILE ( .NOT. LBEND(IDX) ) 1372 WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') 'CC_SETB11>', 1373 & (IBTRAN(I,IDX),I=1,3),(IBDOTS(I,IDX),I=1,MXVEC) 1374 IDX = IDX + 1 1375 END DO 1376 CALL FLSHFO(LUPRI) 1377 CALL QUIT('Overflow error in CC_SETB11.') 1378 END IF 1379 1380 RETURN 1381 END 1382 1383*---------------------------------------------------------------------* 1384* END OF SUBROUTINE CC_SETB11 * 1385*---------------------------------------------------------------------* 1386c /* deck CC_SETD1111 */ 1387*=====================================================================* 1388 SUBROUTINE CC_SETD1111(IDTRAN,IDDOTS,MXTRAN,MXVEC,IZETAV, 1389 & ITAMPA,ITAMPB,ITAMPC,ITAMPD,ITRAN,IVEC) 1390*---------------------------------------------------------------------* 1391* 1392* Purpose: set up list of D matrix transformations 1393* assumes that T^A, T^B, T^C, T^D members of the same list 1394* 1395* IDTRAN - list of D matrix transformations 1396* IDDOTS - list of vectors it should be dottet on 1397* 1398* MXTRAN - maximum list dimension 1399* MXVEC - maximum second dimension for IDDOTS 1400* 1401* IZETAV - index of lagrangian multiplier vector 1402* ITAMPA - index of amplitude vector A 1403* ITAMPB - index of amplitude vector B 1404* ITAMPC - index of amplitude vector C 1405* ITAMPD - index of amplitude vector D 1406* 1407* ITRAN - index in IDTRAN list 1408* IVEC - second index in IDDOTS list 1409* 1410* Written by Christof Haettig, june 1997. 1411* 1412*=====================================================================* 1413 IMPLICIT NONE 1414 1415#include "priunit.h" 1416 INTEGER MXVEC, MXTRAN 1417 INTEGER IDTRAN(5,MXTRAN) 1418 INTEGER IDDOTS(MXVEC,MXTRAN) 1419 1420 LOGICAL LFND 1421 INTEGER IZETAV, ITAMPA, ITAMPB, ITAMPC, ITAMPD 1422 INTEGER ITRAN, IVEC 1423 INTEGER ITAMP, I, IDX 1424 1425* statement functions: 1426 LOGICAL LDTST, LDXTST, LDEND 1427 INTEGER IA, IB, IC, ID, IA1, IB1, IC1, ID1 1428 1429 LDXTST(ITRAN,IA1,IB1,IC1,ID1) = 1430 & IDTRAN(1,ITRAN).EQ.IA1 .AND. IDTRAN(2,ITRAN).EQ.IB1 .AND. 1431 & IDTRAN(3,ITRAN).EQ.IC1 .AND. IDTRAN(4,ITRAN).EQ.ID1 1432 1433 LDTST(ITRAN,IA,IB,IC,ID) = 1434 & LDXTST(ITRAN,IA,IB,IC,ID) .OR. LDXTST(ITRAN,IA,IB,ID,IC) .OR. 1435 & LDXTST(ITRAN,IB,IA,IC,ID) .OR. LDXTST(ITRAN,IB,IA,ID,IC) .OR. 1436 & LDXTST(ITRAN,IA,IC,IB,ID) .OR. LDXTST(ITRAN,IA,IC,ID,IB) .OR. 1437 & LDXTST(ITRAN,IC,IA,IB,ID) .OR. LDXTST(ITRAN,IC,IA,ID,IB) .OR. 1438 & LDXTST(ITRAN,IA,ID,IC,IB) .OR. LDXTST(ITRAN,IA,ID,IB,IC) .OR. 1439 & LDXTST(ITRAN,ID,IA,IC,IB) .OR. LDXTST(ITRAN,ID,IA,IB,IC) .OR. 1440 & LDXTST(ITRAN,ID,IB,IC,IA) .OR. LDXTST(ITRAN,ID,IB,IA,IC) .OR. 1441 & LDXTST(ITRAN,IB,ID,IC,IA) .OR. LDXTST(ITRAN,IB,ID,IA,IC) .OR. 1442 & LDXTST(ITRAN,ID,IC,IB,IA) .OR. LDXTST(ITRAN,ID,IC,IA,IB) .OR. 1443 & LDXTST(ITRAN,IC,ID,IB,IA) .OR. LDXTST(ITRAN,IC,ID,IA,IB) .OR. 1444 & LDXTST(ITRAN,ID,IA,IC,IB) .OR. LDXTST(ITRAN,ID,IA,IB,IC) .OR. 1445 & LDXTST(ITRAN,IA,ID,IC,IB) .OR. LDXTST(ITRAN,IA,ID,IB,IC) 1446 1447 LDEND(ITRAN) = ITRAN.GT.MXTRAN .OR. 1448 & (IDTRAN(1,ITRAN)+IDTRAN(2,ITRAN)+ 1449 & IDTRAN(3,ITRAN)+IDTRAN(4,ITRAN)).LE.0 1450 1451*---------------------------------------------------------------------* 1452* set up list of D matrix transformations 1453*---------------------------------------------------------------------* 1454 ITRAN = 1 1455 LFND = LDTST(ITRAN,ITAMPA,ITAMPB,ITAMPC,ITAMPD) 1456 1457 DO WHILE ( .NOT. (LFND.OR.LDEND(ITRAN)) ) 1458 ITRAN = ITRAN + 1 1459 LFND = LDTST(ITRAN,ITAMPA,ITAMPB,ITAMPC,ITAMPD) 1460 END DO 1461 1462 IF (.NOT.LFND) THEN 1463 IDTRAN(1,ITRAN) = ITAMPA 1464 IDTRAN(2,ITRAN) = ITAMPB 1465 IDTRAN(3,ITRAN) = ITAMPC 1466 IDTRAN(4,ITRAN) = ITAMPD 1467 IDTRAN(5,ITRAN) = 0 1468 ITAMP = IZETAV 1469 ELSE 1470 ITAMP = IZETAV 1471 END IF 1472 1473 IVEC = 1 1474 DO WHILE (IDDOTS(IVEC,ITRAN).NE.ITAMP .AND. 1475 & IDDOTS(IVEC,ITRAN).NE.0 .AND. IVEC.LE.MXVEC) 1476 IVEC = IVEC + 1 1477 END DO 1478 1479 IDDOTS(IVEC,ITRAN) = ITAMP 1480 1481*---------------------------------------------------------------------* 1482 IF (IVEC.GT.MXVEC .OR. ITRAN.GT.MXTRAN) THEN 1483 WRITE (LUPRI,*) 'IVEC :',IVEC 1484 WRITE (LUPRI,*) 'ITRAN:',ITRAN 1485 WRITE (LUPRI,*) 'ITAMPA,ITAMPB,ITAMPC,ITAMPD:', 1486 & ITAMPA,ITAMPB,ITAMPC,ITAMPD 1487 IDX = 1 1488 DO WHILE ( .NOT. LDEND(IDX) ) 1489 WRITE(LUPRI,'(A,5I5,5X,(12I5,20X))') 'CC_SETD1111>', 1490 & (IDTRAN(I,IDX),I=1,5),(IDDOTS(I,IDX),I=1,MXVEC) 1491 IDX = IDX + 1 1492 END DO 1493 CALL FLSHFO(LUPRI) 1494 CALL QUIT('Overflow error in CC_SETD1111.') 1495 END IF 1496 1497 RETURN 1498 END 1499 1500*---------------------------------------------------------------------* 1501* END OF SUBROUTINE CC_SETD1111 * 1502*---------------------------------------------------------------------* 1503c /* deck CC_SETH2112 */ 1504*=====================================================================* 1505 SUBROUTINE CC_SETH2112(IHTRAN,IHDOTS,MXTRAN,MXVEC,IZETAV, 1506 & ITAMPA,ITAMPB,ITAMPC,ITAMPD,ITRAN,IVEC) 1507*---------------------------------------------------------------------* 1508* 1509* Purpose: maintain a list of dot products of H matrix 1510* transformations with right amplitude vectors: 1511* (Z*H*T^A*T^B*T^C) * T^D 1512* N.B.: assumes that T^A and T^D vectors are members of one list 1513* and that T^B and T^C vectors members of a second list 1514* 1515* IHTRAN - list of H matrix transformations 1516* IHDOTS - list of vectors it should be dottet on 1517* 1518* MXTRAN - maximum list dimension 1519* MXVEC - maximum second dimesion for IHDOTS 1520* 1521* IZETAV - index of lagrangian multiplier vector 1522* ITAMPA - index of amplitude vector A 1523* ITAMPB - index of amplitude vector B 1524* ITAMPC - index of amplitude vector C 1525* ITAMPD - index of amplitude vector D 1526* 1527* ITRAN - index in IHTRAN list 1528* IVEC - second index in IHDOTS list 1529* 1530* Written by Christof Haettig, june 1997. 1531* 1532*=====================================================================* 1533 IMPLICIT NONE 1534 1535#include "priunit.h" 1536 INTEGER MXVEC, MXTRAN 1537 INTEGER IHTRAN(5,MXTRAN) 1538 INTEGER IHDOTS(MXVEC,MXTRAN) 1539 1540 LOGICAL LFNDA, LFNDD 1541 INTEGER IZETAV, ITAMPA, ITAMPB, ITAMPC, ITAMPD 1542 INTEGER ITRAN, IVEC 1543 INTEGER ITAMP, I, IDX 1544 1545* statement functions: 1546 LOGICAL LHTST, LHEND 1547 INTEGER IL, IA, IB,IC 1548 LHTST(ITRAN,IL,IA,IB,IC) = IHTRAN(1,ITRAN).EQ.IL .AND. ( 1549 & (IHTRAN(2,ITRAN).EQ.IA .AND. IHTRAN(3,ITRAN).EQ.IB 1550 & .AND. IHTRAN(4,ITRAN).EQ.IC) .OR. 1551 & (IHTRAN(2,ITRAN).EQ.IA .AND. IHTRAN(3,ITRAN).EQ.IC 1552 & .AND. IHTRAN(4,ITRAN).EQ.IB) ) 1553 LHEND(ITRAN) = ITRAN.GT.MXTRAN .OR. 1554 & (IHTRAN(1,ITRAN)+IHTRAN(2,ITRAN)+ 1555 & IHTRAN(3,ITRAN)+IHTRAN(4,ITRAN) ).LE.0 1556 1557*---------------------------------------------------------------------* 1558* set up list of H matrix transformations 1559*---------------------------------------------------------------------* 1560 ITRAN = 1 1561 LFNDA = LHTST(ITRAN,IZETAV,ITAMPD,ITAMPB,ITAMPC) 1562 LFNDD = LHTST(ITRAN,IZETAV,ITAMPA,ITAMPB,ITAMPC) 1563 1564 DO WHILE ( .NOT. (LFNDA.OR.LFNDD.OR.LHEND(ITRAN)) ) 1565 ITRAN = ITRAN + 1 1566 LFNDA = LHTST(ITRAN,IZETAV,ITAMPD,ITAMPB,ITAMPC) 1567 LFNDD = LHTST(ITRAN,IZETAV,ITAMPA,ITAMPB,ITAMPC) 1568 END DO 1569 1570 IF (.NOT.(LFNDA.OR.LFNDD)) THEN 1571 IHTRAN(1,ITRAN) = IZETAV 1572 IHTRAN(2,ITRAN) = ITAMPA 1573 IHTRAN(3,ITRAN) = ITAMPB 1574 IHTRAN(4,ITRAN) = ITAMPC 1575 IHTRAN(5,ITRAN) = 0 1576 ITAMP = ITAMPD 1577 ELSE 1578 IF (LFNDA) ITAMP = ITAMPA 1579 IF (LFNDD) ITAMP = ITAMPD 1580 END IF 1581 1582 IVEC = 1 1583 DO WHILE (IHDOTS(IVEC,ITRAN).NE.ITAMP .AND. 1584 & IHDOTS(IVEC,ITRAN).NE.0 .AND. IVEC.LE.MXVEC) 1585 IVEC = IVEC + 1 1586 END DO 1587 1588 IHDOTS(IVEC,ITRAN) = ITAMP 1589 1590*---------------------------------------------------------------------* 1591 IF (IVEC.GT.MXVEC .OR. ITRAN.GT.MXTRAN) THEN 1592 WRITE (LUPRI,*) 'IVEC :',IVEC 1593 WRITE (LUPRI,*) 'ITRAN:',ITRAN 1594 WRITE (LUPRI,*) 'IZETAV,ITAMPA,ITAMPB,ITAMPC,ITAMPD:', 1595 & IZETAV,ITAMPA,ITAMPB,ITAMPC,ITAMPD 1596 IDX = 1 1597 DO WHILE (.NOT. LHEND(IDX)) 1598 WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') 'CC_SETH2112>', 1599 & (IHTRAN(I,IDX),I=1,4),(IHDOTS(I,IDX),I=1,MXVEC) 1600 IDX = IDX + 1 1601 END DO 1602 CALL FLSHFO(LUPRI) 1603 CALL QUIT('Overflow error in CC_SETH2112.') 1604 END IF 1605 1606 RETURN 1607 END 1608*---------------------------------------------------------------------* 1609* END OF SUBROUTINE CC_SETH2112 * 1610*---------------------------------------------------------------------* 1611c /* deck CC_SETC211 */ 1612*=====================================================================* 1613 SUBROUTINE CC_SETC211(ICTRAN,ICDOTS,MXTRAN,MXVEC, 1614 & IZETAV,ITAMPA,ITAMPB,ITAMPC,ITRAN,IVEC) 1615*---------------------------------------------------------------------* 1616* 1617* Purpose: set up list of C matrix transformations 1618* assumes that T^A vectors on one list and 1619* T^B and T^C vectors members of a second list 1620* 1621* ICTRAN - list of C matrix transformations 1622* ICDOTS - list of vectors it should be dottet on 1623* 1624* MXTRAN - maximum list dimension 1625* MXVEC - maximum second dimension for ICDOTS 1626* 1627* IZETAV - index of lagrangian multiplier vector 1628* ITAMPA - index of amplitude vector A 1629* ITAMPB - index of amplitude vector B 1630* ITAMPC - index of amplitude vector C 1631* 1632* ITRAN - index in ICTRAN list 1633* IVEC - second index in ICDOTS list 1634* 1635* Written by Christof Haettig, june 1997. 1636* 1637*=====================================================================* 1638 IMPLICIT NONE 1639 1640#include "priunit.h" 1641 INTEGER MXVEC, MXTRAN 1642 INTEGER ICTRAN(4,MXTRAN) 1643 INTEGER ICDOTS(MXVEC,MXTRAN) 1644 1645 LOGICAL LFND 1646 INTEGER IZETAV, ITAMPA, ITAMPB, ITAMPC 1647 INTEGER ITRAN, IVEC 1648 INTEGER ITAMP, I, IDX 1649 1650* statement functions: 1651 LOGICAL LCTST, LCEND 1652 INTEGER IB, IA, IC 1653 LCTST(ITRAN,IA,IB,IC) = ICTRAN(1,ITRAN).EQ.IA .AND. ( 1654 & ( ICTRAN(2,ITRAN).EQ.IB .AND. ICTRAN(3,ITRAN).EQ.IC ) .OR. 1655 & ( ICTRAN(2,ITRAN).EQ.IC .AND. ICTRAN(3,ITRAN).EQ.IB ) ) 1656 LCEND(ITRAN) = ITRAN.GT.MXTRAN .OR. 1657 & (ICTRAN(1,ITRAN)+ICTRAN(2,ITRAN)+ICTRAN(3,ITRAN)).LE.0 1658 1659*---------------------------------------------------------------------* 1660* set up list of B matrix transformations 1661*---------------------------------------------------------------------* 1662 ITRAN = 1 1663 LFND = LCTST(ITRAN,ITAMPA,ITAMPB,ITAMPC) 1664 1665 DO WHILE ( .NOT. (LFND.OR.LCEND(ITRAN)) ) 1666 ITRAN = ITRAN + 1 1667 LFND = LCTST(ITRAN,ITAMPA,ITAMPB,ITAMPC) 1668 END DO 1669 1670 IF (.NOT.LFND) THEN 1671 ICTRAN(1,ITRAN) = ITAMPA 1672 ICTRAN(2,ITRAN) = ITAMPB 1673 ICTRAN(3,ITRAN) = ITAMPC 1674 ICTRAN(4,ITRAN) = 0 1675 ITAMP = IZETAV 1676 ELSE 1677 ITAMP = IZETAV 1678 END IF 1679 1680 IVEC = 1 1681 DO WHILE (ICDOTS(IVEC,ITRAN).NE.ITAMP .AND. 1682 & ICDOTS(IVEC,ITRAN).NE.0 .AND. IVEC.LE.MXVEC) 1683 IVEC = IVEC + 1 1684 END DO 1685 1686 IF (.FALSE. .AND. ICDOTS(IVEC,ITRAN).NE.ITAMP) THEN 1687 WRITE (LUPRI,'(A,5I5)') 'CC_SETC211> add ', 1688 & ITAMPA,ITAMPB,ITAMPC,IZETAV 1689 END IF 1690 1691 ICDOTS(IVEC,ITRAN) = ITAMP 1692 1693*---------------------------------------------------------------------* 1694 IF (IVEC.GT.MXVEC .OR. ITRAN.GT.MXTRAN) THEN 1695 WRITE (LUPRI,*) 'IVEC :',IVEC 1696 WRITE (LUPRI,*) 'ITRAN:',ITRAN 1697 WRITE (LUPRI,*) 'ITAMPA,ITAMPB,ITAMPC:',ITAMPA,ITAMPB,ITAMPC 1698 IDX = 1 1699 DO WHILE ( .NOT. LCEND(IDX) ) 1700 WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') 'CC_SETC211>', 1701 & (ICTRAN(I,IDX),I=1,4),(ICDOTS(I,IDX),I=1,MXVEC) 1702 IDX = IDX + 1 1703 END DO 1704 CALL FLSHFO(LUPRI) 1705 CALL QUIT('Overflow error in CC_SETC211.') 1706 END IF 1707 1708 RETURN 1709 END 1710 1711*---------------------------------------------------------------------* 1712* END OF SUBROUTINE CC_SETC211 * 1713*---------------------------------------------------------------------* 1714