diff options
Diffstat (limited to 'Dragon/src/SHIDRV.f')
| -rw-r--r-- | Dragon/src/SHIDRV.f | 443 |
1 files changed, 443 insertions, 0 deletions
diff --git a/Dragon/src/SHIDRV.f b/Dragon/src/SHIDRV.f new file mode 100644 index 0000000..be16d07 --- /dev/null +++ b/Dragon/src/SHIDRV.f @@ -0,0 +1,443 @@ +*DECK SHIDRV + SUBROUTINE SHIDRV (IPLIB,IPTRK,IFTRAK,LEVEL,NGRO,NBISO,NBMIX, + 1 NREG,NUN,CDOOR,NRES,IMPX,ISONRF,ISONAM,MIX,DEN,SN,SB,LSHI, + 2 IPHASE,IPROB,MAT,VOL,KEYFLX,LEAKSW,TITR,IGRMIN,IGRMAX,MAXX0, + 3 IBIEFF,IGC,ITRANZ,EPS) +* +*----------------------------------------------------------------------- +* +*Purpose: +* Perform a multidimensional self-shielding calculation in order to +* compute the dilution cross section of each resonant isotope present +* in the domain. +* +*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 +* IPLIB pointer to the internal microscopic cross section library +* (L_LIBRARY signature). +* IPTRK pointer to the tracking. (L_TRACK signature). +* IFTRAK unit number of the sequential binary tracking file. +* LEVEL type of self-shielding model (=0 original Stamm'ler model; +* =1 original Stamm'ler model with Nordheim approximation; +* =2 Stamm'ler model with Nordheim approximation and Riemann +* integration method). +* NGRO number of energy groups. +* NBISO number of isotopes present in the calculation domain. +* NBMIX number of mixtures in the macrolib. +* NREG number of regions. +* NUN number of unknowns in the flux or source vector in one +* energy group. +* CDOOR name of the geometry/solution module. +* NRES number of resonant mixtures. +* IMPX print flag. +* ISONRF reference name of isotopes. +* ISONAM alias name of isotopes. +* MIX mix number of each isotope (can be zero). +* DEN density of each isotope. +* LSHI resonant region number associated with each isotope. +* Infinite dilution will be assumed if LSHI(i)=0. +* IPHASE type of flux solution (=1 use a native flux solution door; +* =2 use collision probabilities). +* IPROB adjoint macrolib flag (=0 direct; =1 adjoint). +* MAT index-number of the mixture type assigned to each volume. +* VOL volumes. +* KEYFLX pointers of fluxes in unknown vector. +* LEAKSW leakage flag (=.TRUE. only if leakage is present on the outer +* surface). +* TITR title. +* IGRMIN first group where the self-shielding is applied. +* IGRMAX most thermal group where the self-shielding is applied. +* MAXX0 maximum number of self-shielding iterations. +* IBIEFF Livolant-Jeanpierre normalization flag (=1 to activate). +* IGC Goldstein-Cohen approximation flag (=1 to activate). +* The Goldstein-Cohen approximation is activated only in cases +* where they are available in the internal library. +* ITRANZ type of transport correction used in the self-shielding +* calculations. +* EPS convergence criterion for the self-shielding iterations. +* +*Parameters: input/output +* SN on input, estimate of the dilution cross section in each +* energy group of each isotope. A value of 1.0e10 is used +* for infinite dilution and; +* on output, computed dilution cross section in each energy +* group of each isotope. +* +*Parameters: output +* SB dilution cross section as used in Livolant-Jeanpierre +* normalization. +* +*----------------------------------------------------------------------- +* + USE GANLIB +*---- +* SUBROUTINE ARGUMENTS +*---- + TYPE(C_PTR) IPLIB,IPTRK + INTEGER IFTRAK,LEVEL,NGRO,NBISO,NBMIX,NREG,NUN,NRES,IMPX, + 1 ISONRF(3,NBISO),ISONAM(3,NBISO),MIX(NBISO),LSHI(NBISO), + 2 IPHASE,IPROB,MAT(NREG),KEYFLX(NREG),IGRMIN,IGRMAX,MAXX0,IBIEFF, + 3 IGC,ITRANZ + REAL DEN(NBISO),SN(NGRO,NBISO),SB(NGRO,NBISO),VOL(NREG),EPS + LOGICAL LEAKSW + CHARACTER CDOOR*12,TITR*72 +*---- +* LOCAL VARIABLES +*---- + PARAMETER (NALPHA=9,NRAT=(NALPHA+1)/2,NSTATE=40) + TYPE(C_PTR) JPLIB,KPLIB + INTEGER IPAR(NSTATE) + REAL TMPDAY(3) + CHARACTER HSMG*130 + LOGICAL START,BIEFF,LGC,LOGDO +*---- +* ALLOCATABLE ARRAYS +*---- + INTEGER, ALLOCATABLE, DIMENSION(:) :: LSHI2 + REAL, ALLOCATABLE, DIMENSION(:) :: SIGE + REAL, ALLOCATABLE, DIMENSION(:,:) :: SIGT1,SIGT2,SIGT3 + LOGICAL, ALLOCATABLE, DIMENSION(:) :: MASK,MASKL + LOGICAL, ALLOCATABLE, DIMENSION(:,:) :: NOCONV +*---- +* SCRATCH STORAGE ALLOCATION +*---- + ALLOCATE(SIGT1(NBMIX,NGRO),SIGT2(NBMIX,NGRO),SIGT3(NBMIX,NGRO)) + ALLOCATE(MASK(NBMIX),MASKL(NGRO),NOCONV(NBMIX,NGRO)) +* + IF(IMPX.GE.5) THEN + WRITE (6,'(//23H SHIDRV: VALUES OF MAT:)') + I1=1 + KI=(NREG-1)/11+1 + DO 10 I=1,KI + I2=I1+10 + IF(I2.GT.NREG) I2=NREG + WRITE (6,340) (J,J=I1,I2) + WRITE (6,350) (MAT(J),J=I1,I2) + I1=I1+11 + 10 CONTINUE + WRITE (6,'(//)') + ENDIF +*---- +* RECOVER SELF SHIELDING DATA +*---- + IF(LEAKSW) CALL XABORT('SHIDRV: NEUTRON LEAKAGE IS FORBIDDEN.') + IF(CDOOR.EQ.' ') CALL XABORT('SHIDRV: THE GEOMETRY IS NOT YET ' + 1 //'DEFINED.') + BIEFF=(IBIEFF.EQ.1) + LGC=(IGC.EQ.1) +* + IF(IMPX.GT.0) THEN + WRITE (6,400) TITR,CDOOR + WRITE (6,'(25H STAMM''LER APPROXIMATION./)') + WRITE (6,405) IGRMIN,IGRMAX,MAXX0,IBIEFF,IGC,ITRANZ,LEVEL, + 1 IPHASE,EPS + ENDIF + IF(NRES.EQ.0) THEN + WRITE (6,410) + RETURN + ENDIF + DO 30 I=1,NREG + IF(MAT(I).GT.NBMIX) THEN + WRITE (HSMG,380) NBMIX + CALL XABORT(HSMG) + ENDIF + 30 CONTINUE + IGRMAX=MIN(IGRMAX,NGRO) + DO 55 LLL=1,NGRO + DO 40 I=1,NBISO + SB(LLL,I)=SN(LLL,I) + 40 CONTINUE + DO 50 IBM=1,NBMIX + NOCONV(IBM,LLL)=.FALSE. + 50 CONTINUE + 55 CONTINUE + CALL LCMSIX(IPLIB,'MACROLIB',1) + JPLIB=LCMGID(IPLIB,'GROUP') + DO 70 LLL=IGRMIN,IGRMAX + IF(IPROB.EQ.0) LL=LLL + IF(IPROB.EQ.1) LL=NGRO-LLL+1 + KPLIB=LCMGIL(JPLIB,LL) + CALL LCMGET(KPLIB,'NTOT0',SIGT2(1,LLL)) +*---- +* TRANSPORT CORRECTION +*---- + IF(ITRANZ.NE.0) THEN + CALL LCMGET(KPLIB,'TRANC',SIGT3(1,LLL)) + ELSE + SIGT3(:NBMIX,LLL)=0.0 + ENDIF +* + NOCONV(:NBMIX,LLL)=.TRUE. + 70 CONTINUE + CALL LCMSIX(IPLIB,' ',2) + IF(IMPX.GE.5) THEN + WRITE (6,'(/20H SHIBA INPUT VALUES:/)') + DO 80 L=IGRMIN,IGRMAX + WRITE(6,420) L + WRITE(6,460) (SN(L,J),J=1,NBISO) + WRITE(6,480) (SIGT2(IBM,L),IBM=1,NBMIX) + 80 CONTINUE + WRITE(6,490) + ENDIF +*---- +* ELIMINATE ISOTOPE ABSENT FROM GEOMETRY +*---- + DO IBM=1,NBMIX + DO IREG=1,NREG + IF(MAT(IREG).EQ.IBM) GO TO 85 + ENDDO + DO ISO=1,NBISO + IF(MIX(ISO).EQ.IBM) LSHI(ISO)=0 + ENDDO + 85 CONTINUE + ENDDO +*---- +* RECOMPUTE THE VECTOR LSHI +*---- + IF(LEVEL.EQ.0) THEN + NRES2=NRES + ELSE + ALLOCATE(LSHI2(NBISO)) + NRES1=0 + NRES2=0 + DO 90 ISO=1,NBISO + LSHI2(ISO)=0 + 90 CONTINUE + DO 140 INRS=1,NRES + 100 DENMAX=0.0 + KSOT=0 + DO 120 ISO=1,NBISO + IF(LSHI2(ISO).EQ.0) THEN + VOLISO=0.0 + DO 110 I=1,NREG + IF(MAT(I).EQ.MIX(ISO)) VOLISO=VOLISO+VOL(I) + 110 CONTINUE + IF((LSHI(ISO).EQ.INRS).AND.(DEN(ISO)*VOLISO.GT.DENMAX)) THEN + KSOT=ISO + DENMAX=DEN(ISO)*VOLISO + ENDIF + ENDIF + 120 CONTINUE + IF(KSOT.GT.0) THEN + NRES2=NRES2+1 + DO 130 ISO=1,NBISO + IF((ISONRF(1,ISO).EQ.ISONRF(1,KSOT)).AND. + 1 (ISONRF(2,ISO).EQ.ISONRF(2,KSOT)).AND. + 2 (ISONRF(3,ISO).EQ.ISONRF(3,KSOT)).AND. + 3 (LSHI(ISO).EQ.INRS)) LSHI2(ISO)=NRES2 + IF((ISONAM(1,ISO).EQ.ISONAM(1,KSOT)).AND. + 1 (ISONAM(2,ISO).EQ.ISONAM(2,KSOT)).AND. + 2 (LSHI(ISO).EQ.INRS)) LSHI2(ISO)=NRES2 + 130 CONTINUE + GO TO 100 + ENDIF + IF(NRES2.EQ.NRES1) THEN + WRITE(HSMG,'(43HSHIDRV: NO RESONANT ISOTOPES IN RESONANT RE, + 1 11HGION NUMBER,I4,5H (1).)') INRS + CALL XABORT(HSMG) + ENDIF + NRES1=NRES2 + 140 CONTINUE + ENDIF +*---- +* DETERMINE THE AMOUNT OF SCRATCH STORAGE REQUIRED +*---- + NBMIX2=0 + DO 150 ISO=1,NBISO + IF(LSHI(ISO).GT.0) NBMIX2=NBMIX2+1 + 150 CONTINUE +*---- +* ITERATION LOOP +*---- + IF(LEVEL.EQ.0) ALLOCATE(SIGE(NRES2*NGRO)) + NITER=0 + 160 NITER=NITER+1 + START=(NITER.EQ.1) + IF(IMPX.GT.5) WRITE (6,430) NITER + DO 175 L=IGRMIN,IGRMAX + DO 170 IBM=1,NBMIX + SIGT1(IBM,L)=SIGT2(IBM,L) + 170 CONTINUE + 175 CONTINUE + IF(LEVEL.EQ.0) THEN + CALL SHISN2 (IPLIB,IPTRK,IFTRAK,NGRO,NBISO,NBMIX,NREG,NUN, + 1 CDOOR,NRES,NBMIX2,IMPX,ISONAM,MIX,DEN,SN,SB,LSHI,IPHASE, + 2 MAT,VOL,KEYFLX,LEAKSW,TITR,START,SIGT2,SIGT3,NOCONV,BIEFF, + 3 LGC,SIGE) + ELSE + DO 210 INRS=1,NRES2 + NBNRS=0 + DO 200 IBM=1,NBMIX + LOGDO=.FALSE. + DO 180 I=1,NREG + LOGDO=LOGDO.OR.(MAT(I).EQ.IBM) + 180 CONTINUE + IF(.NOT.LOGDO) GO TO 200 + DO 190 ISO=1,NBISO + IF((MIX(ISO).EQ.IBM).AND.(LSHI2(ISO).EQ.INRS)) THEN + NBNRS=NBNRS+1 + GO TO 200 + ENDIF + 190 CONTINUE + 200 CONTINUE + IF(NBNRS.EQ.0) THEN + IF(START.AND.(IMPX.GE.1)) WRITE(6,385) 'SHIDRV',INRS + GO TO 210 + ELSE IF(START.AND.(NBNRS.GT.1).AND.(IMPX.GE.1)) THEN + WRITE (6,370) NBNRS,INRS + ENDIF + CALL SHISN3 (IPLIB,IPTRK,IFTRAK,LEVEL,NGRO,NBISO,NBMIX,NREG, + 1 NUN,CDOOR,INRS,NBNRS,IMPX,ISONAM,MIX,DEN,SN,SB,LSHI2,IPHASE, + 2 MAT,VOL,KEYFLX,LEAKSW,TITR,START,SIGT2,SIGT3,NOCONV,BIEFF,LGC) + 210 CONTINUE + ENDIF + ZZMAX=0.0 + LNGRO=0 + ICOUNT=0 + DO 240 L=IGRMIN,IGRMAX + ZNORM=0.0 + DO 220 IBM=1,NBMIX + ZNORM=MAX(ZNORM,ABS(SIGT2(IBM,L))) + 220 CONTINUE + ZMAX=0.0 + MASKL(L)=.FALSE. + DO 230 IBM=1,NBMIX + YMAX=ABS(SIGT1(IBM,L)-SIGT2(IBM,L))/ZNORM + ZMAX=MAX(ZMAX,YMAX) + NOCONV(IBM,L)=(NOCONV(IBM,L).AND.(YMAX.GT.EPS)) + MASKL(L)=MASKL(L).OR.NOCONV(IBM,L) + 230 CONTINUE + IF(MASKL(L)) ICOUNT=ICOUNT+1 + IF(ZMAX.GT.ZZMAX) THEN + ZZMAX=ZMAX + LNGRO=L + ENDIF + 240 CONTINUE + IF(IMPX.GE.2) WRITE (6,440) NITER,ICOUNT,ZZMAX,LNGRO + IF(IMPX.GE.10) THEN + WRITE (6,450) (L,MASKL(L),L=IGRMIN,IGRMAX) + WRITE (6,'(/31H INPUT MACROSCOPIC X-S IN GROUP,I4,1H:)') LNGRO + WRITE (6,'(1X,1P,10E12.4)') (SIGT1(IBM,LNGRO),IBM=1,NBMIX) + WRITE (6,'(/32H OUTPUT MACROSCOPIC X-S IN GROUP,I4,1H:)') LNGRO + WRITE (6,'(1X,1P,10E12.4)') (SIGT2(IBM,LNGRO),IBM=1,NBMIX) + ENDIF + IF(IMPX.GT.3) THEN + WRITE (6,'(/29H OUTPUT DILUTION X-S IN GROUP,I4,1H:)') LNGRO + WRITE (6,'(1X,1P,10E12.4)') (SN(LNGRO,J),J=1,NBISO) + ENDIF + IF((NITER.GE.MAXX0).AND.(ICOUNT.GT.0)) THEN + WRITE (6,390) + GO TO 250 + ELSE IF(ICOUNT.GT.0) THEN + GO TO 160 + ENDIF +*---- +* CONVERGENCE IS OBTAINED +*---- + 250 IF(LEVEL.GT.0) DEALLOCATE(LSHI2) + IF(LEVEL.EQ.0) DEALLOCATE(SIGE) + IF(IMPX.GE.3) THEN + WRITE (6,'(/21H SHIBA OUTPUT VALUES:/)') + DO 260 L=IGRMIN,IGRMAX + WRITE(6,420) L + WRITE(6,460) (SN(L,J),J=1,NBISO) + IF(BIEFF) WRITE(6,470) (SB(L,J),J=1,NBISO) + IF(IMPX.GE.5) WRITE(6,480) (SIGT2(IBM,L),IBM=1,NBMIX) + 260 CONTINUE + WRITE(6,490) + ENDIF +*---- +* COMPUTE THE NEW SELF-SHIELDED MACROSCOPIC CROSS SECTIONS +*---- + MASKL(:NGRO)=.FALSE. + DO 280 LLL=IGRMIN,IGRMAX + MASKL(LLL)=.TRUE. + 280 CONTINUE + DO 300 IBM=1,NBMIX + DO 290 ISO=1,NBISO + MASK(IBM)=(MIX(ISO).EQ.IBM).AND.(LSHI(ISO).GT.0) + IF(MASK(IBM)) GO TO 300 + 290 CONTINUE + 300 CONTINUE + ITSTMP=0 + TMPDAY(1)=0.0 + TMPDAY(2)=0.0 + TMPDAY(3)=0.0 + CALL LIBMIX(IPLIB,NBMIX,NGRO,NBISO,ISONAM,MIX,DEN,MASK,MASKL, + > ITSTMP,TMPDAY) + IF(IMPX.GT.0) WRITE (6,500) NITER,ZZMAX +*---- +* STORE THE GENERAL SHIBA PARAMETERS +*---- + IPAR(:NSTATE)=0 + IPAR(1)=IGRMIN + IPAR(2)=IGRMAX + IPAR(3)=MAXX0 + IPAR(4)=IBIEFF + IPAR(5)=IGC + IPAR(6)=ITRANZ + IPAR(7)=LEVEL + IPAR(8)=IPHASE + CALL LCMSIX(IPLIB,'SHIBA',1) + CALL LCMPUT(IPLIB,'STATE-VECTOR',NSTATE,1,IPAR) + CALL LCMPUT(IPLIB,'EPS-SHIBA',1,2,EPS) + CALL LCMSIX(IPLIB,' ',2) +*---- +* SCRATCH STORAGE DEALLOCATION +*---- + DEALLOCATE(NOCONV,MASKL,MASK) + DEALLOCATE(SIGT3,SIGT2,SIGT1) + RETURN +* + 340 FORMAT(//26H VOLUME NB. :,11(I5,3X,1HI)) + 350 FORMAT( 26H MIXTURE (MAT) :,11(I5,3X,1HI)) + 370 FORMAT(/42H SHIDRV: USE THE NORDHEIM MODEL TO PROCESS,I3,5H RESO, + 1 39HNANT MIXTURES IN RESONANT REGION NUMBER,I3,1H.) + 380 FORMAT(32HSHIDRV: INVALID VALUE OF NBMIX (,I5,2H).) + 385 FORMAT(A6,1X,': RESONANT REGION =',I10,1X,'NOT USED.') + 390 FORMAT(/1X,61(1H*)/42H SHIDRV: MAXIMUM NUMBER OF SELF-SHIELDING , + 1 20HITERATIONS EXCEEDED./1X,61(1H*)/) + 400 FORMAT( + > 1X,'SHIBA MULTIDIMENSIONAL SELF-SHIELDING CALCULATION', + > 1X,'-> A. HEBERT'/ + > 1X,A72/ + > 1X,'COLLISION PROBABILITY MODULE: ',A12/) + 405 FORMAT(/8H OPTIONS/8H -------/ + 1 7H IGRMIN,I8,27H (FIRST GROUP TO PROCESS)/ + 2 7H IGRMAX,I8,34H (MOST THERMAL GROUP TO PROCESS)/ + 3 7H MAXX0 ,I8,33H (MAXIMUM NUMBER OF ITERATIONS)/ + 4 7H IBIEFF,I8,46H (=1: USE LIVOLANT-JEANPIERRE NORMALIZATION)/ + 5 7H IGC ,I8,42H (=1: USE GOLDSTEIN-COHEN APPROXIMATION)/ + 6 7H ITRANZ,I8,45H (0=NO TRANSPORT CORRECTION/1=APOLLO TYPE/2, + 7 57H=RECOVER FROM LIBRARY/3=WIMS-D TYPE/4=LEAKAGE CORRECTION)/ + 8 7H LEVEL ,I8,46H (=0: STAMM'LER; =1: STAMM'LER/NORDHEIM; =2:, + 9 18H RIEMANN/NORDHEIM)/ + 1 7H IPHASE,I8,37H (=1: NATIVE ASSEMBLY; =2: USE PIJ)/ + 2 7H EPS ,1P,E8.1,22H (STOPING CRITERION)/) + 410 FORMAT(/52H SHIDRV: THERE IS NO REQUEST TO PROCESS ANY RESONANT, + 1 9H ISOTOPE./) + 420 FORMAT(/1X,131(1H-)//8H GROUP =,I4/) + 430 FORMAT(/40H PERFORMING SELF-SHIELDING ITERATION NB.,I5) + 440 FORMAT(/27H SELF-SHIELDING ITERATION =,I4,5X,14HNUMBER OF NON , + 1 18HCONVERGED GROUPS =,I4,5X,7HERROR =,1P,E13.4,0P,9H IN GROUP, + 2 I4/) + 450 FORMAT(7H MASKL(,I3,2H)=,L1,:,8H MASKL(,I3,2H)=,L1,:,8H MASKL(, + 1 I3,2H)=,L1,:,8H MASKL(,I3,2H)=,L1,:,8H MASKL(,I3,2H)=,L1,:, + 2 8H MASKL(,I3,2H)=,L1,:,8H MASKL(,I3,2H)=,L1,:,8H MASKL(,I3, + 3 2H)=,L1,:,8H MASKL(,I3,2H)=,L1) + 460 FORMAT(/37H MICROSCOPIC DILUTION CROSS SECTIONS:/(9X,1P,11E11.3)) + 470 FORMAT(/53H LIVOLANT AND JEANPIERRE MICROSCOPIC DILUTION CROSS S, + 1 8HECTIONS:/(9X,1P,11E11.3)) + 480 FORMAT(/34H MACROSCOPIC TOTAL CROSS SECTIONS:/(9X,1P,11E11.3)) + 490 FORMAT(/1X,131(1H-)/) + 500 FORMAT(/41H CONVERGENCE REACHED AT SHIBA ITERATION =,I4,7H ERROR, + 1 2H =,1P,E11.3/) + END |