From 7dfcc480ba1e19bd3232349fc733caef94034292 Mon Sep 17 00:00:00 2001 From: stainer_t Date: Mon, 8 Sep 2025 13:48:49 +0200 Subject: Initial commit from Polytechnique Montreal --- Dragon/src/LIBDRB.f | 773 ++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 773 insertions(+) create mode 100644 Dragon/src/LIBDRB.f (limited to 'Dragon/src/LIBDRB.f') diff --git a/Dragon/src/LIBDRB.f b/Dragon/src/LIBDRB.f new file mode 100644 index 0000000..600f737 --- /dev/null +++ b/Dragon/src/LIBDRB.f @@ -0,0 +1,773 @@ +*DECK LIBDRB + SUBROUTINE LIBDRB (IPDRL,NGRO,NL,NDEL,NBESP,SN,SB,NED,HVECT,DELTA, + 1 LBIN,NFS,BENER,AWR,DELECC,IGECCO,IMPX,NGF,NGFR,LSCAT,LSIGF,LADD, + 2 LGOLD,SIGS,SCAT,TOTAL,ZNPHI,SIGF,CHI,CHI4G,SADD,GOLD,BIN) +* +*----------------------------------------------------------------------- +* +*Purpose: +* Read an interpolate in dilution one isotope in draglib format at a +* selected temperature. +* +*Copyright: +* Copyright (C) 2002 Ecole Polytechnique de Montreal +* This library is free software; you can redistribute it and/or +* modify it under the terms of the GNU Lesser General Public +* License as published by the Free Software Foundation; either +* version 2.1 of the License, or (at your option) any later version +* +*Author(s): A. Hebert +* +*Parameters: input +* IPDRL pointer to the draglib (L_DRAGLIB signature). +* NGRO number of energy groups. +* NL number of Legendre orders required in the calculation. +* NL=1 or higher. +* NDEL number of delayed precursor groups. +* NBESP number of energy-dependent fission spectra. +* SN dilution cross section in each energy group. A value of +* 1.0E10 is used for infinite dilution. +* SB dilution cross section as used in Livolant and Jeanpierre +* normalization. +* NED number of extra vector edits. +* HVECT names of the extra vector edits. +* DELTA lethargy widths. +* LBIN number of fine groups. +* NFS number of fine groups per coarse group. +* BENER energy limits of the fine groups. +* AWR mass ratio for current isotope. +* DELECC lethargy width of eccolib libraries. +* IGECCO number of equal-width lethargy groups with eccolib libraries. +* IMPX print flag. +* +*Parameters: input/output +* NGF number of fast groups without self-shielding. +* NGFR number of fast and resonance groups. +* +*Parameters: output +* LSCAT scattering mask (=.true. if a given Legendre order of the +* scattering cross section exists). +* LSIGF fission mask (=.true. if the isotope can fission). +* LADD additional cross section mask (=.true. if a given additional +* cross section exists). +* LGOLD Goldstein-Cohen mask (=.true. if Goldstein-Cohen parameters +* exists). +* SIGS scattering cross sections. +* SCAT scattering transfer matrices. +* TOTAL total cross sections. +* ZNPHI fluxes. +* SIGF nu*fission cross sections. +* CHI fission spectrum. +* CHI4G energy-dependent fission spectra. +* SADD additional cross sections. +* GOLD Goldstein-Cohen parameters. +* BIN BIN(IGR,1): total fine group cross sections; +* BIN(IGR,2): isotropic scattering fine group cross sections; +* BIN(IGR,3): nu*fission fine group cross sections. +* +*----------------------------------------------------------------------- +* + USE GANLIB +*---- +* SUBROUTINE ARGUMENTS +*---- + PARAMETER(MAXDIL=50) + CHARACTER*(*) HVECT(NED) + TYPE(C_PTR) IPDRL + INTEGER NGRO,NL,NDEL,NBESP,NED,LBIN,NFS(NGRO),IGECCO,IMPX,NGF,NGFR + REAL SN(NGRO),SB(NGRO),DELTA(NGRO),BENER(LBIN+1),AWR,DELECC, + 1 SIGS(NGRO,NL),SCAT(NGRO,NGRO,NL),TOTAL(NGRO),ZNPHI(NGRO), + 2 SIGF(NGRO,0:NDEL),CHI(NGRO,0:NDEL),CHI4G(NGRO,NBESP), + 3 SADD(NGRO,NED),GOLD(NGRO),BIN(LBIN,3) + LOGICAL LSCAT(NL),LSIGF,LADD(NED),LGOLD +*---- +* LOCAL VARIABLES +*---- + CHARACTER CM*2,CD*4,HSMG*131,HNUSIG*12,HCHI*12,HTOTAL*5 + PARAMETER (IOUT=6,MAXTRA=10000) + INTEGER KTOTLR,KSIGFR,KCHIR,KPHIR + LOGICAL LPCAT + DOUBLE PRECISION TMP,ZNGAR,SQD,SQ0,SQ1,SQ2,SQ3,FACT1,FACT2 + REAL DILUT(MAXDIL) +*---- +* ALLOCATABLE ARRAYS +*---- + INTEGER, ALLOCATABLE, DIMENSION(:) :: NJJ,IJJ,KADDR + REAL, ALLOCATABLE, DIMENSION(:) :: GAR,PRI,STIS,UUU,SSS + REAL, ALLOCATABLE, DIMENSION(:,:) :: TERP,SIGT + LOGICAL, ALLOCATABLE, DIMENSION(:) :: LSDIL,LPDIL,LINF +*---- +* SCRATCH STORAGE ALLOCATION +*---- + ALLOCATE(NJJ(NGRO),IJJ(NGRO),KADDR(NED)) + ALLOCATE(GAR(NGRO**2),TERP(MAXDIL,NGRO),SIGT(NGRO,MAXDIL)) + ALLOCATE(LSDIL(NL),LPDIL(NL),LINF(NGRO)) +* + TOTAL(:NGRO)=0.0 + HTOTAL='NTOT0' + CALL LCMLEN(IPDRL,'NTOT0',LENGT,ITYLCM) + IF(LENGT.EQ.0) CALL XABORT('LIBDRB: MISSING TOTAL XS INFO.') + CALL LCMGET(IPDRL,HTOTAL,TOTAL) + CALL LCMLEN(IPDRL,'NUSIGF',LENGT,ITYLCM) + LSIGF=(LENGT.GT.0) + IF(LSIGF) THEN + DO 10 IDEL=0,NDEL + IF(IDEL.EQ.0) THEN + HNUSIG='NUSIGF' + HCHI='CHI' + ELSE + WRITE(HNUSIG,'(6HNUSIGF,I2.2)') IDEL + WRITE(HCHI,'(3HCHI,I2.2)') IDEL + CALL LCMLEN(IPDRL,HNUSIG,ILONG,ITYLCM) + IF(ILONG.EQ.0) CALL XABORT('LIBDRB: MISSING '//HNUSIG// + 1 ' INFO.') + ENDIF + SIGF(:NGRO,IDEL)=0.0 + CALL LCMGET(IPDRL,HNUSIG,SIGF(1,IDEL)) + IF((NBESP.EQ.0).OR.(IDEL.GT.0)) THEN + CHI(:NGRO,IDEL)=0.0 + CALL LCMLEN(IPDRL,HCHI,LENGT,ITYLCM) + IF(LENGT.GT.0) CALL LCMGET(IPDRL,HCHI,CHI(1,IDEL)) + ENDIF + 10 CONTINUE + DO 15 ISP=1,NBESP + WRITE(HCHI,'(5HCHI--,I2.2)') ISP + CHI4G(:NGRO,ISP)=0.0 + CALL LCMLEN(IPDRL,HCHI,LENGT,ITYLCM) + IF(LENGT.GT.0) CALL LCMGET(IPDRL,HCHI,CHI4G(1,ISP)) + 15 CONTINUE + ELSE + SIGF(:NGRO,0:NDEL)=0.0 + ENDIF + DO 150 IL=0,NL-1 + SIGS(:NGRO,IL+1)=0.0 + SCAT(:NGRO,:NGRO,IL+1)=0.0 + WRITE (CM,'(I2.2)') IL + CALL LCMLEN(IPDRL,'SCAT'//CM,LENGT,ITYLCM) + LPCAT=(LENGT.GT.0) + IF(LPCAT.AND.(IGECCO.EQ.0)) THEN + CALL LCMGET(IPDRL,'NJJS'//CM,NJJ) + CALL LCMGET(IPDRL,'IJJS'//CM,IJJ) + LENGT2=0 + DO 20 I=1,NGRO + LENGT2=LENGT2+NJJ(I) + 20 CONTINUE + GAR(:LENGT)=0.0 + CALL LCMGET(IPDRL,'SCAT'//CM,GAR) + IGAR=0 +* IG2 IS THE SECONDARY GROUP. + DO 40 IG2=1,NGRO + DO 30 IG1=IJJ(IG2),IJJ(IG2)-NJJ(IG2)+1,-1 + IGAR=IGAR+1 + SCAT(IG2,IG1,IL+1)=GAR(IGAR) + 30 CONTINUE + 40 CONTINUE + ELSE IF(LPCAT) THEN + ! on-flight elastic scattering kernel + CALL LCMGET(IPDRL,'NJJS'//CM,NJJ) + CALL LCMGET(IPDRL,'IJJS'//CM,IJJ) + ALLOCATE(PRI(MAXTRA),STIS(NGRO),UUU(NGRO),SSS(IGECCO)) + CALL LIBPRI(MAXTRA,DELECC,AWR,0,IL,NPRI,PRI) + LENGT2=0 + DO 50 I=1,NGRO + LENGT2=LENGT2+NJJ(I) + 50 CONTINUE + GAR(:LENGT)=0.0 + CALL LCMGET(IPDRL,'SCAT'//CM,GAR) + UUU(1)=DELTA(1) + DO 60 I=2,NGRO + UUU(I)=UUU(I-1)+DELTA(I) + 60 CONTINUE + IGAR=0 +* IG2 IS THE SECONDARY GROUP. + DO 90 IG2=1,NGRO + IF(IG2.LE.IGECCO) THEN + CALL LIBECT(MAXTRA,IG2,PRI,UUU,DELECC,DELTA,NPRI,1,MML,STIS) + IGAR=IGAR+NJJ(IG2) + SSS(IG2)=GAR(IGAR) + DO 70 I=1,MML + IG1=IG2-I+1 + IF(IG1.LE.0) GO TO 90 + SCAT(IG2,IG1,IL+1)=STIS(I)*SSS(IG1) + 70 CONTINUE + ELSE + DO 80 IG1=IJJ(IG2),IJJ(IG2)-NJJ(IG2)+1,-1 + IGAR=IGAR+1 + SCAT(IG2,IG1,IL+1)=GAR(IGAR) + 80 CONTINUE + ENDIF + 90 CONTINUE + DEALLOCATE(SSS,UUU,STIS,PRI) + ENDIF + CALL LCMLEN(IPDRL,'SIGS'//CM,LENGT,ITYLCM) + LSCAT(IL+1)=(LENGT.GT.0) + IF(LSCAT(IL+1)) THEN + CALL LCMGET(IPDRL,'SIGS'//CM,SIGS(1,IL+1)) + CALL LCMLEN(IPDRL,'PCAT'//CM,LENGT,ITYLCM) + IF(.NOT.LPCAT.AND.(LENGT.GT.0)) THEN + CALL LCMGET(IPDRL,'NJJS'//CM,NJJ) + CALL LCMGET(IPDRL,'IJJS'//CM,IJJ) + LENGT=0 + DO 100 I=1,NGRO + LENGT=LENGT+NJJ(I) + 100 CONTINUE + GAR(:LENGT)=0.0 + CALL LCMGET(IPDRL,'PCAT'//CM,GAR) + IGAR=0 + DO 115 IG2=1,NGRO + DO 110 IG1=IJJ(IG2),IJJ(IG2)-NJJ(IG2)+1,-1 + IGAR=IGAR+1 + SCAT(IG2,IG1,IL+1)=GAR(IGAR)*SIGS(IG1,IL+1) + 110 CONTINUE + 115 CONTINUE + ELSE IF(.NOT.LPCAT) THEN + DO 120 IG1=1,NGRO + SCAT(IG1,IG1,IL+1)=SIGS(IG1,IL+1) + 120 CONTINUE + ENDIF + ELSE IF(LPCAT) THEN + DO 140 IG1=1,NGRO + TMP=0.0D0 + DO 130 IG2=1,NGRO + TMP=TMP+SCAT(IG2,IG1,IL+1) + 130 CONTINUE + SIGS(IG1,IL+1)=REAL(TMP) + 140 CONTINUE + LSCAT(IL+1)=.TRUE. + ENDIF + 150 CONTINUE + LSCAT(1)=.TRUE. + DO 160 IED=1,NED + SADD(:NGRO,IED)=0.0 + CALL LCMLEN(IPDRL,HVECT(IED),LENGT,ITYLCM) + LADD(IED)=(LENGT.GT.0) + IF(LADD(IED)) CALL LCMGET(IPDRL,HVECT(IED),SADD(1,IED)) + 160 CONTINUE + CALL LCMLEN(IPDRL,'NGOLD',LENGT,ITYLCM) + LGOLD=(LENGT.GT.0) + IF(LGOLD) THEN + GOLD(:NGRO)=0.0 + CALL LCMGET(IPDRL,'NGOLD',GOLD) + ELSE + GOLD(:NGRO)=1.0 + ENDIF + IF(LBIN.GT.0) THEN + CALL LCMGET(IPDRL,'BIN-'//HTOTAL,BIN(1,1)) + CALL LCMGET(IPDRL,'BIN-SIGS00',BIN(1,2)) + CALL LCMLEN(IPDRL,'BIN-NUSIGF',LENGF,ITYLCM) + IF(LENGF.GT.0) THEN + CALL LCMGET(IPDRL,'BIN-NUSIGF',BIN(1,3)) + ELSE + BIN(:LBIN,3)=0.0 + ENDIF + IGF0=0 + DO 190 IG=1,NGRO + IF(NFS(IG).GT.0) THEN +* BIN CROSS SECTION NORMALIZATION. + SQ0=0.0D0 + SQ1=0.0D0 + SQ2=0.0D0 + SQ3=0.0D0 + DO 170 IGF=IGF0+1,IGF0+NFS(IG) + DELTAU=LOG(BENER(IGF)/BENER(IGF+1)) + SQ0=SQ0+DELTAU + SQ1=SQ1+BIN(IGF,1)*DELTAU + SQ2=SQ2+(BIN(IGF,1)-BIN(IGF,2))*DELTAU + SQ3=SQ3+BIN(IGF,3)*DELTAU + 170 CONTINUE + FACT1=TOTAL(IG)*(SQ0/SQ1) + FACT2=(TOTAL(IG)-SIGS(IG,1))*(SQ0/SQ2) + DO 180 IGF=IGF0+1,IGF0+NFS(IG) + BIN(IGF,2)=REAL(BIN(IGF,2)*FACT2+BIN(IGF,1)*(FACT1-FACT2)) + BIN(IGF,1)=REAL(BIN(IGF,1)*FACT1) + IF((LENGF.GT.0).AND.(SQ3.NE.0.0)) THEN + BIN(IGF,3)=REAL(BIN(IGF,3)*(SIGF(IG,0)*(SQ0/SQ3))) + ENDIF + 180 CONTINUE + IGF0=IGF0+NFS(IG) + ENDIF + 190 CONTINUE + ENDIF + KTOTLR=0 + KSIGFR=0 + KPHIR=0 + KCHIR=0 + KADR=0 + KADDR(:NED)=0 +*---- +* PERFORM DILUTION INTERPOLATION. +*---- + CALL LCMLEN(IPDRL,'DILUTION',NDIL,ITYLCM) + IF(NDIL.GT.0) THEN + IF(NDIL+1.GT.MAXDIL) CALL XABORT('LIBDRB: INVALID MAXDIL.') + CALL LCMGET(IPDRL,'DILUTION',DILUT) + IF(DILUT(NDIL).GE.1.0E10) CALL XABORT('LIBDRB: INVALID DILUTI' + 1 //'ON VALUE.') +*---- +* FIND MAX LENGTH OF VECTORS ON SUBMAT +* KTOTLR,KSIGFR,KCHIR,KPHIR AND KADDR +* GIVES LENGTH OF SELF SHIELDING VECTOR +* FOR TOTAL, SIGF, CHI, PHI AND ADD XS +*---- + LSDIL(:NL)=.FALSE. + LPDIL(:NL)=.FALSE. + DO 240 IDIL=1,NDIL + WRITE (CD,'(I4.4)') IDIL + CALL LCMSIX(IPDRL,'SUBMAT'//CD,1) + SIGT(:NGRO,IDIL)=0.0 + CALL LCMGET(IPDRL,HTOTAL,SIGT(1,IDIL)) + DO 220 IL=0,NL-1 + WRITE (CM,'(I2.2)') IL + CALL LCMLEN(IPDRL,'SCAT'//CM,LENGT,ITYLCM) + IF(.NOT.LSDIL(IL+1)) + > LSDIL(IL+1)=(LENGT.GT.0).AND.LSCAT(IL+1) + CALL LCMLEN(IPDRL,'SIGS'//CM,LENGT,ITYLCM) + IF(.NOT.LPDIL(IL+1)) + > LPDIL(IL+1)=(LENGT.GT.0).AND.LSCAT(IL+1) + 220 CONTINUE + CALL LCMLEN(IPDRL,HTOTAL,LENGT,ITYLCM) + KTOTLR=MAX(KTOTLR,LENGT) + CALL LCMLEN(IPDRL,'NUSIGF',LENGT,ITYLCM) + KSIGFR=MAX(KSIGFR,LENGT) + IF(NBESP.EQ.0) THEN + CALL LCMLEN(IPDRL,'CHI',LENGT,ITYLCM) + KCHIR=MAX(KCHIR,LENGT) + ELSE + DO 225 ISP=1,NBESP + WRITE(HCHI,'(5HCHI--,I2.2)') ISP + CALL LCMLEN(IPDRL,HCHI,LENGT,ITYLCM) + KCHIR=MAX(KCHIR,LENGT) + 225 CONTINUE + ENDIF + CALL LCMLEN(IPDRL,'NWT0',LENGT,ITYLCM) + KPHIR=MAX(KPHIR,LENGT) + DO 230 IED=1,NED + CALL LCMLEN(IPDRL,HVECT(IED),LENGT,ITYLCM) + IF((LENGT.GT.0).AND.LADD(IED)) THEN + KADDR(IED)=MAX(KADDR(IED),LENGT) + KADR=MAX(KADDR(IED),KADR) + ENDIF + 230 CONTINUE + CALL LCMSIX(IPDRL,' ',2) + 240 CONTINUE + NGRRE=MAX(KTOTLR,KSIGFR,KCHIR,KPHIR,KADR) + IF(NGRRE.GT.NGRO) CALL XABORT('LIBDRB: TOO MANY GROUPS.') +* + TERP(:MAXDIL,:NGRO)=0.0 + DILUT(NDIL+1)=1.0E10 + DO 280 IG1=1,NGRRE + LINF(IG1)=.FALSE. + ZNPHI(IG1)=0.0 + DILX=MIN(SN(IG1),1.0E10) + IF(DILX.LE.0.0) THEN + WRITE (HSMG,930) IG1 + CALL XABORT(HSMG) + ENDIF + IFIRST=0 + DO 260 I=1,NDIL+1 + IF(ABS(DILX-DILUT(I)).LE.1.0E-5*ABS(DILX)) THEN + TERP(I,IG1)=1.0 + GO TO 280 + ELSE IF(DILX.LT.DILUT(I)) THEN + IFIRST=I-1 + GO TO 270 + ENDIF + 260 CONTINUE +* + 270 SQD=SQRT(DILX) + IF((IFIRST-1.GE.1).AND.(IFIRST+2.LE.NDIL)) THEN + SQ0=SQRT(DILUT(IFIRST-1)) + SQ1=SQRT(DILUT(IFIRST)) + SQ2=SQRT(DILUT(IFIRST+1)) + SQ3=SQRT(DILUT(IFIRST+2)) + TERP(IFIRST-1,IG1)=REAL((SQ1-SQD)*(SQ2-SQD)*(SQ3-SQD)/ + 1 (SQ1-SQ0)/(SQ2-SQ0)/(SQ3-SQ0)) + TERP(IFIRST,IG1)=REAL((SQ0-SQD)*(SQ2-SQD)*(SQ3-SQD)/ + 1 (SQ0-SQ1)/(SQ2-SQ1)/(SQ3-SQ1)) + TERP(IFIRST+1,IG1)=REAL((SQ0-SQD)*(SQ1-SQD)*(SQ3-SQD)/ + 1 (SQ0-SQ2)/(SQ1-SQ2)/(SQ3-SQ2)) + TERP(IFIRST+2,IG1)=REAL((SQ0-SQD)*(SQ1-SQD)*(SQ2-SQD)/ + 1 (SQ0-SQ3)/(SQ1-SQ3)/(SQ2-SQ3)) + TT=TERP(IFIRST-1,IG1)*SIGT(IG1,IFIRST-1) + 1 +TERP(IFIRST,IG1)*SIGT(IG1,IFIRST) + 2 +TERP(IFIRST+1,IG1)*SIGT(IG1,IFIRST+1) + 3 +TERP(IFIRST+2,IG1)*SIGT(IG1,IFIRST+2) + YMIN=MIN(SIGT(IG1,IFIRST),SIGT(IG1,IFIRST+1)) + YMAX=MAX(SIGT(IG1,IFIRST),SIGT(IG1,IFIRST+1)) + IF((TT.GT.YMAX).OR.(TT.LT.YMIN)) THEN + TERP(IFIRST-1,IG1)=0.0 + TERP(IFIRST,IG1)=REAL((SQ2-SQD)/(SQ2-SQ1)) + TERP(IFIRST+1,IG1)=REAL((SQ1-SQD)/(SQ1-SQ2)) + TERP(IFIRST+2,IG1)=0.0 + ENDIF + ELSE IF((IFIRST.EQ.1).AND.(IFIRST+2.LE.NDIL)) THEN + SQ1=SQRT(DILUT(1)) + SQ2=SQRT(DILUT(2)) + SQ3=SQRT(DILUT(3)) + TERP(1,IG1)=REAL((SQ2-SQD)*(SQ3-SQD)/(SQ2-SQ1)/(SQ3-SQ1)) + TERP(2,IG1)=REAL((SQ1-SQD)*(SQ3-SQD)/(SQ1-SQ2)/(SQ3-SQ2)) + TERP(3,IG1)=REAL((SQ1-SQD)*(SQ2-SQD)/(SQ1-SQ3)/(SQ2-SQ3)) + TT=TERP(1,IG1)*SIGT(IG1,1)+TERP(2,IG1)*SIGT(IG1,2) + 1 +TERP(3,IG1)*SIGT(IG1,3) + YMIN=MIN(SIGT(IG1,1),SIGT(IG1,2)) + YMAX=MAX(SIGT(IG1,1),SIGT(IG1,2)) + IF((TT.GT.YMAX).OR.(TT.LT.YMIN)) THEN + TERP(1,IG1)=REAL((SQ2-SQD)/(SQ2-SQ1)) + TERP(2,IG1)=REAL((SQ1-SQD)/(SQ1-SQ2)) + TERP(3,IG1)=0.0 + ENDIF + ELSE IF((IFIRST-1.GE.1).AND.(IFIRST+1.EQ.NDIL)) THEN + SQ0=SQRT(DILUT(NDIL-2)) + SQ1=SQRT(DILUT(NDIL-1)) + SQ2=SQRT(DILUT(NDIL)) + TERP(NDIL-2,IG1)=REAL((SQ1-SQD)*(SQ2-SQD)/(SQ1-SQ0)/(SQ2-SQ0)) + TERP(NDIL-1,IG1)=REAL((SQ0-SQD)*(SQ2-SQD)/(SQ0-SQ1)/(SQ2-SQ1)) + TERP(NDIL,IG1)=REAL((SQ0-SQD)*(SQ1-SQD)/(SQ0-SQ2)/(SQ1-SQ2)) + TT=TERP(NDIL-2,IG1)*SIGT(IG1,NDIL-2) + 1 +TERP(NDIL-1,IG1)*SIGT(IG1,NDIL-1) + 2 +TERP(NDIL,IG1)*SIGT(IG1,NDIL) + YMIN=MIN(SIGT(IG1,NDIL-1),SIGT(IG1,NDIL)) + YMAX=MAX(SIGT(IG1,NDIL-1),SIGT(IG1,NDIL)) + IF((TT.GT.YMAX).OR.(TT.LT.YMIN)) THEN + TERP(NDIL-2,IG1)=0.0 + TERP(NDIL-1,IG1)=REAL((SQ2-SQD)/(SQ2-SQ1)) + TERP(NDIL,IG1)=REAL((SQ1-SQD)/(SQ1-SQ2)) + ENDIF + ELSE IF((IFIRST.EQ.0).OR.((IFIRST.EQ.1).AND.(NDIL.EQ.2))) THEN + SQ0=SQRT(DILUT(1)) + SQ1=SQRT(DILUT(2)) + TERP(1,IG1)=REAL((SQ1-SQD)/(SQ1-SQ0)) + TERP(2,IG1)=REAL((SQ0-SQD)/(SQ0-SQ1)) + ELSE IF(IFIRST.EQ.NDIL) THEN + LINF(IG1)=.TRUE. + TERP(NDIL,IG1)=DILUT(NDIL)/DILX + ELSE + CALL XABORT('LIBDRB: FAILURE OF DILUTION INTERPOLATION.') + ENDIF + 280 CONTINUE +* + NGRODP=NGRO+1 + NGROIN=0 + DO 330 IDIL=1,NDIL + NCORF=0 + DO 290 IG1=NGRO,1,-1 + IF(TERP(IDIL,IG1).NE.0.0) THEN + NCORF=IG1 + GO TO 300 + ENDIF + 290 CONTINUE + 300 NGROIN=MAX(NCORF,NGROIN) + NCORD=NGRO+1 + DO 310 IG1=1,NGROIN + IF(TERP(IDIL,IG1).NE.0.0) THEN + NCORD=IG1 + GO TO 320 + ENDIF + 310 CONTINUE + 320 NGRODP=MIN(NCORD,NGRODP) + 330 CONTINUE + DO 345 IDIL=1,NDIL + DO 340 IG1=1,NGRO + IF(SIGT(IG1,IDIL).NE.0.0) THEN + NGF=MIN(NGF,IG1-1) + NGFR=MAX(NGFR,IG1) + ENDIF + 340 CONTINUE + 345 CONTINUE + IF(NGROIN.EQ.0.OR.NGRRE.EQ.0) THEN + ZNPHI(:NGRO)=1.0 + GO TO 850 + ENDIF + KTOTLR=MIN(KTOTLR,NGROIN) + KSIGFR=MIN(KSIGFR,NGROIN) + KCHIR=MIN(KCHIR,NGROIN) + KPHIR=MIN(KPHIR,NGROIN) + DO 360 IED=1,NED + KADDR(IED)=MIN(KADDR(IED),NGROIN) + 360 CONTINUE +*---- +* VARIOUS DIMENSION OF VECTORS ARE SET +* LOOP OVER DILUTION AND SELF-SHIELD XS +* FROM NGRODP TO NGROIN (THESE CORRESPOND +* TO CASES WHERE DIL<1.0E10 FOR AT LEAST ONE GROUP +* HERE ONE ASSUMES THAT TOTAL XS ALWAYS SELF-SHIELDED +*---- + DO 640 IDIL=1,NDIL + DO 370 IG1=1,NGRO + IF(TERP(IDIL,IG1).NE.0.0) GO TO 380 + 370 CONTINUE + GO TO 640 + 380 WRITE (CD,'(I4.4)') IDIL + CALL LCMSIX(IPDRL,'SUBMAT'//CD,1) + DO 390 IG1=NGRODP,NGROIN + TOTAL(IG1)=TOTAL(IG1)+TERP(IDIL,IG1)*SIGT(IG1,IDIL) + 390 CONTINUE + IF(KSIGFR.GT.0) THEN + DO 420 IDEL=0,NDEL + IF(IDEL.EQ.0) THEN + HNUSIG='NUSIGF' + ELSE + WRITE(HNUSIG,'(6HNUSIGF,I2.2)') IDEL + ENDIF + GAR(:KSIGFR)=0.0 + CALL LCMGET(IPDRL,HNUSIG,GAR) + DO 410 IG1=NGRODP,KSIGFR + SIGF(IG1,IDEL)=SIGF(IG1,IDEL)+TERP(IDIL,IG1)*GAR(IG1) + 410 CONTINUE + 420 CONTINUE + ENDIF + IF(KCHIR.GT.0) THEN + DO 440 IDEL=0,NDEL + IF(IDEL.EQ.0) THEN + IF(NBESP.GT.0) GO TO 440 + HCHI='CHI' + ELSE + WRITE(HCHI,'(3HCHI,I2.2)') IDEL + ENDIF + GAR(:KCHIR)=0.0 + CALL LCMGET(IPDRL,HCHI,GAR) + DO 430 IG1=NGRODP,KCHIR + CHI(IG1,IDEL)=CHI(IG1,IDEL)+TERP(IDIL,IG1)*GAR(IG1) + 430 CONTINUE + 440 CONTINUE + DO 460 ISP=1,NBESP + WRITE(HCHI,'(5HCHI--,I2.2)') ISP + GAR(:KCHIR)=0.0 + CALL LCMLEN(IPDRL,HCHI,ILONG,ITYLCM) + IF(ILONG.GT.0) THEN + CALL LCMGET(IPDRL,HCHI,GAR) + DO 450 IG1=NGRODP,KCHIR + CHI4G(IG1,ISP)=CHI4G(IG1,ISP)+TERP(IDIL,IG1)*GAR(IG1) + 450 CONTINUE + ENDIF + 460 CONTINUE + ENDIF + DO 600 IL=0,NL-1 + WRITE (CM,'(I2.2)') IL + IF(LSDIL(IL+1).AND.(IGECCO.EQ.0)) THEN + CALL LCMGET(IPDRL,'NJJS'//CM,NJJ) + CALL LCMGET(IPDRL,'IJJS'//CM,IJJ) + LENGT=0 + DO 470 I=1,NGRO + LENGT=LENGT+NJJ(I) + 470 CONTINUE + GAR(:LENGT)=0.0 + CALL LCMGET(IPDRL,'SCAT'//CM,GAR) + IGAR=0 + DO 490 IG2=1,NGRO + DO 480 IG1=IJJ(IG2),IJJ(IG2)-NJJ(IG2)+1,-1 + IGAR=IGAR+1 + SCAT(IG2,IG1,IL+1)=SCAT(IG2,IG1,IL+1)+TERP(IDIL,IG1) + 1 *GAR(IGAR) + 480 CONTINUE + 490 CONTINUE + ELSE IF(LSDIL(IL+1)) THEN + ! on-flight elastic scattering kernel + CALL LCMGET(IPDRL,'NJJS'//CM,NJJ) + CALL LCMGET(IPDRL,'IJJS'//CM,IJJ) + ALLOCATE(PRI(MAXTRA),STIS(NGRO),UUU(NGRO),SSS(IGECCO)) + CALL LIBPRI(MAXTRA,DELECC,AWR,0,IL,NPRI,PRI) + LENGT2=0 + DO 500 I=1,NGRO + LENGT2=LENGT2+NJJ(I) + 500 CONTINUE + GAR(:LENGT)=0.0 + CALL LCMGET(IPDRL,'SCAT'//CM,GAR) + UUU(1)=DELTA(1) + DO 510 I=2,NGRO + UUU(I)=UUU(I-1)+DELTA(I) + 510 CONTINUE + IGAR=0 +* IG2 IS THE SECONDARY GROUP. + DO 540 IG2=1,NGRO + IF(IG2.LE.IGECCO) THEN + CALL LIBECT(MAXTRA,IG2,PRI,UUU,DELECC,DELTA,NPRI,1,MML, + 1 STIS) + IGAR=IGAR+NJJ(IG2) + SSS(IG2)=GAR(IGAR) + DO 520 I=1,MML + IG1=IG2-I+1 + IF(IG1.LE.0) GO TO 540 + SCAT(IG2,IG1,IL+1)=SCAT(IG2,IG1,IL+1)+TERP(IDIL,IG1)* + 1 STIS(I)*SSS(IG1) + 520 CONTINUE + ELSE + DO 530 IG1=IJJ(IG2),IJJ(IG2)-NJJ(IG2)+1,-1 + IGAR=IGAR+1 + SCAT(IG2,IG1,IL+1)=SCAT(IG2,IG1,IL+1)+TERP(IDIL,IG1)* + 1 GAR(IGAR) + 530 CONTINUE + ENDIF + 540 CONTINUE + DEALLOCATE(SSS,UUU,STIS,PRI) + ENDIF + IF(LPDIL(IL+1)) THEN + GAR(:NGRO)=0.0 + CALL LCMGET(IPDRL,'SIGS'//CM,GAR) + DO 590 IG1=NGRODP,NGROIN + SIGS(IG1,IL+1)=SIGS(IG1,IL+1)+TERP(IDIL,IG1)*GAR(IG1) + 590 CONTINUE + ENDIF + 600 CONTINUE + IF(KPHIR.GT.0) THEN + GAR(:KPHIR)=0.0 + CALL LCMGET(IPDRL,'NWT0',GAR) + DO 610 IG1=NGRODP,KPHIR + IF(.NOT.LINF(IG1)) THEN + ZNPHI(IG1)=ZNPHI(IG1)+TERP(IDIL,IG1)*GAR(IG1)* + 1 DILUT(IDIL) + ELSE + ZNPHI(IG1)=GAR(IG1)*DILUT(IDIL) + ENDIF + 610 CONTINUE + ENDIF + DO 630 IED=1,NED + IF(KADDR(IED).GT.0) THEN + GAR(:KADDR(IED))=0.0 + CALL LCMGET(IPDRL,HVECT(IED),GAR) + DO 620 IG1=NGRODP,KADDR(IED) + SADD(IG1,IED)=SADD(IG1,IED)+TERP(IDIL,IG1)*GAR(IG1) + 620 CONTINUE + ENDIF + 630 CONTINUE + CALL LCMSIX(IPDRL,' ',2) + 640 CONTINUE +*---- +* COMPUTE MISSING SCATTERING INFORMATION. +*---- + DO 710 IL=0,NL-1 + IF(LPDIL(IL+1).AND.(.NOT.LSDIL(IL+1))) THEN + WRITE (CM,'(I2.2)') IL + CALL LCMGET(IPDRL,'NJJS'//CM,NJJ) + CALL LCMGET(IPDRL,'IJJS'//CM,IJJ) + LENGT=0 + DO 650 I=1,NGRO + LENGT=LENGT+NJJ(I) + 650 CONTINUE + GAR(:LENGT)=0.0 + CALL LCMGET(IPDRL,'PCAT'//CM,GAR) + IGAR=0 + DO 680 IG2=1,NGRO + DO 660 IG1=1,NGRO + SCAT(IG2,IG1,IL+1)=0.0 + 660 CONTINUE + DO 670 IG1=IJJ(IG2),IJJ(IG2)-NJJ(IG2)+1,-1 + IGAR=IGAR+1 + SCAT(IG2,IG1,IL+1)=GAR(IGAR)*SIGS(IG1,IL+1) + 670 CONTINUE + 680 CONTINUE + ELSE IF((.NOT.LPDIL(IL+1)).AND.LSDIL(IL+1)) THEN + DO 700 IG1=1,NGRO + TMP=0.0D0 + DO 690 IG2=1,NGRO + TMP=TMP+SCAT(IG2,IG1,IL+1) + 690 CONTINUE + IF(IL.EQ.0) THEN + SIGS(IG1,1)=MIN(REAL(TMP),TOTAL(IG1)) + ELSE + SIGS(IG1,IL+1)=REAL(TMP) + ENDIF + 700 CONTINUE + ENDIF + 710 CONTINUE +*---- +* COMPUTE CONDENSED FINE STRUCTURE FUNCTION. +*---- + DO 740 IG1=1,NGROIN + IF((.NOT.LSDIL(1)).AND.(.NOT.LPDIL(1))) THEN +* SCATTERING CROSS SECTIONS ARE NOT SELF-SHIELDED. + TMP=-TOTAL(IG1) + DO 720 IG2=1,IG1-1 + TMP=TMP+SCAT(IG1,IG2,1)*ZNPHI(IG2)*DELTA(IG2)/DELTA(IG1) + 720 CONTINUE + ZNGAR=(TMP+SCAT(IG1,IG1,1))*SB(IG1)/ + 1 (SB(IG1)-SCAT(IG1,IG1,1)) + ELSE +* SCATTERING CROSS SECTIONS ARE SELF-SHIELDED. + ZNGAR=-TOTAL(IG1) + DO 730 IG2=1,IG1 + ZNGAR=ZNGAR+SCAT(IG1,IG2,1)*DELTA(IG2)/DELTA(IG1) + 730 CONTINUE + ENDIF + IF(IG1.LT.NGRODP) ZNGAR=0.0 + IF(KPHIR.EQ.0) THEN +* USE A CALCULATED VALUE. + ZNPHI(IG1)=REAL(1.0+ZNGAR/SB(IG1)) + ELSE IF(LINF(IG1)) THEN +* USE AN INTERPOLATED VALUE NEAR INFINITE DILUTION. + AUX=(DILUT(NDIL)/SB(IG1))**2 + ZNPHI(IG1)=REAL(AUX*ZNPHI(IG1)+(1.0-AUX)*ZNGAR) + ZNPHI(IG1)=1.0+ZNPHI(IG1)/SB(IG1) + ELSE +* USE AN INTERPOLATED VALUE. + ZNPHI(IG1)=1.0+ZNPHI(IG1)/SB(IG1) + ENDIF + IF((ZNPHI(IG1).LE.0.0).OR.(ZNPHI(IG1).GT.10.0)) THEN + WRITE (HSMG,960) ZNPHI(IG1),IG1,SB(IG1),SN(IG1),KPHIR + CALL XABORT(HSMG) + ELSE IF((ZNPHI(IG1).GT.1.2).AND.(IMPX.GE.3)) THEN + WRITE (HSMG,960) ZNPHI(IG1),IG1,SB(IG1),SN(IG1),KPHIR + WRITE(6,'(1X,A)') HSMG + ENDIF + 740 CONTINUE + DO 750 IG1=NGROIN+1,NGRO + ZNPHI(IG1)=1.0 + 750 CONTINUE +*---- +* DIVIDE EFFECTIVE REACTION RATES BY ZNPHI FOR SELF-SHIELDED +* REACTION RATES +*---- + DO 780 IL=0,NL-1 + IF(LSCAT(IL+1).AND.(LSDIL(IL+1).OR.LPDIL(IL+1))) THEN + DO 770 IG1=NGRODP,NGROIN + SIGS(IG1,IL+1)=SIGS(IG1,IL+1)/ZNPHI(IG1) + DO 760 IG2=1,NGRO + SCAT(IG2,IG1,IL+1)=SCAT(IG2,IG1,IL+1)/ZNPHI(IG1) + 760 CONTINUE + 770 CONTINUE + ENDIF + 780 CONTINUE + DO 790 IG1=NGRODP,NGROIN + TOTAL(IG1)=TOTAL(IG1)/ZNPHI(IG1) + 790 CONTINUE + IF(KSIGFR.GT.0) THEN + DO 810 IDEL=0,NDEL + DO 800 IG1=NGRODP,NGROIN + SIGF(IG1,IDEL)=SIGF(IG1,IDEL)/ZNPHI(IG1) + 800 CONTINUE + 810 CONTINUE + ENDIF + DO 830 IED=1,NED + IF(KADDR(IED).GT.0) THEN + DO 820 IG1=NGRODP,NGROIN + SADD(IG1,IED)=SADD(IG1,IED)/ZNPHI(IG1) + 820 CONTINUE + ENDIF + 830 CONTINUE + IF(IMPX.GT.4) THEN + WRITE(IOUT,940) + DO 840 IG1=1,NGRO + WRITE (IOUT,950) IG1,SN(IG1),SB(IG1),ZNPHI(IG1),TOTAL(IG1), + 1 SIGS(IG1,1),SIGF(IG1,0),GOLD(IG1) + 840 CONTINUE + WRITE (IOUT,'(/)') + ENDIF + ELSE + ZNPHI(:NGRO)=1.0 + ENDIF +*---- +* SCRATCH STORAGE DEALLOCATION +*---- + 850 DEALLOCATE(LINF,LPDIL,LSDIL) + DEALLOCATE(SIGT,TERP,GAR) + DEALLOCATE(KADDR,IJJ,NJJ) + RETURN +* + 930 FORMAT(42HLIBDRB: NEGATIVE OR ZERO DILUTION IN GROUP,I4,1H.) + 940 FORMAT(/5X,'GROUP',10X,'DILUT',13X,'SB',11X,'NWT0',10X,'NTOT0', + 1 11X,'SIGS',9X,'NUSIGF',10X,'NGOLD') + 950 FORMAT(5X,I5,1P,8E15.5) + 960 FORMAT(32HLIBDRB: INVALID VALUE OF ZNPHI (,1P,E11.3, + 1 10H) IN GROUP,I4,11H. DILUTION=,E11.3,2H (,E11.3, + 2 9H). KPHIR=,I4,1H.) + END -- cgit v1.2.3