*DECK DOORPV SUBROUTINE DOORPV (CDOOR,JPSYS,NPSYS,IPTRK,IFTRAK,IMPX,NGRP,NREG, 1 NBMIX,NANI,MAT,VOL,KNORM,IPIJK,LEAKSW,ITPIJ,LNORM,TITR,NALBP) * *----------------------------------------------------------------------- * *Purpose: * Calculation of the collision probabilities. Vectorial version. * *Copyright: * Copyright (C) 2004 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 * CDOOR name of the geometry/solution operator. * JPSYS pointer to the PIJ LCM object (L_PIJ signature). JPSYS is * a list of directories. * NPSYS index array pointing to the JPSYS list component corresponding * to each energy group. Set to zero if a group is not to be * processed. Usually, NPSYS(I)=I. * IPTRK pointer to the tracking (L_TRACK signature). * IFTRAK unit of the sequential binary tracking file. * IMPX print flag (equal to zero for no print). * NGRP number of energy groups. * NREG total number of merged blocks for which specific values * of the neutron flux and reactions rates are required. * NBMIX number of mixtures (NBMIX=max(MAT(i))). * NANI number of Legendre orders (usually equal to one). * MAT index-number of the mixture type assigned to each volume. * VOL volumes. * KNORM normalization scheme for PIJ matrices. * IPIJK pij option (=1 pij, =4 pijk). * LEAKSW leakage flag (=.true. if neutron leakage through external * boundary is present). * ITPIJ type of collision probability information available: * =1 scattering modified pij (wij); * =2 standard pij; * =3 scattering modified pij+pijk (wij,wijk); * =4 standard pij+pijk. * LNORM logical switch for removing leakage from collision * probabilities and keeping the PIS information. * TITR title. * NALBP number of physical albedos. * *----------------------------------------------------------------------- * USE GANLIB *---- * SUBROUTINE ARGUMENTS *---- CHARACTER CDOOR*12,TITR*72 LOGICAL LEAKSW,LNORM TYPE(C_PTR) JPSYS,IPTRK INTEGER NPSYS(NGRP),IFTRAK,IMPX,NGRP,NREG,NBMIX,NANI,MAT(NREG), > KNORM,IPIJK,ITPIJ,NALBP REAL VOL(NREG) INTEGER NNPSYS(1) *---- * LOCAL VARIABLES *---- PARAMETER(NSTATE=40,IOUT=6) CHARACTER CMS(3)*1 INTEGER ISTATE(NSTATE) LOGICAL LBIHET TYPE(C_PTR) KPSYS INTEGER, ALLOCATABLE, DIMENSION(:) :: MAT2 TYPE(C_PTR) :: PREG_PTR *---- * ALLOCATABLE ARRAYS *---- INTEGER, ALLOCATABLE, DIMENSION(:) :: MATALB REAL, ALLOCATABLE, DIMENSION(:) :: VOL2,PGAR,SGAR,SGAS,PROBKS, > PIS REAL, ALLOCATABLE, DIMENSION(:,:) :: PGARG,ALBP REAL, POINTER, DIMENSION(:,:) :: PREG REAL, ALLOCATABLE, DIMENSION(:,:) :: VOLSUR DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: DPROB,DPROBX *---- * DATA STATEMENT AND INLINE FUNCTION *---- SAVE CMS DATA CMS/'1','2','3'/ INDPOS(I,J)=MAX(I,J)*(MAX(I,J)-1)/2+MIN(I,J) *---- * DOUBLE HETEROGENEITY TREATMENT *---- NNPSYS(1)=1 CALL LCMGET(IPTRK,'STATE-VECTOR',ISTATE) IF(ISTATE(1).NE.NREG) CALL XABORT('DOORPV: WRONG VALUE OF NREG.') NSOUT=ISTATE(5) LBIHET=ISTATE(40).NE.0 NREGAR=0 NBMIXG=0 IF(LBIHET) THEN ALLOCATE(MAT2(NREG),VOL2(NREG)) DO I=1,NREG MAT2(I)=MAT(I) VOL2(I)=VOL(I) ENDDO NREGAR=NREG NBMIXG=NBMIX CALL DOORAB(CDOOR,JPSYS,NPSYS,IPTRK,IMPX,NGRP,NREG,NBMIX,NANI, 1 MAT,VOL) ENDIF * NELPIJ=NREG*(NREG+1)/2 NB1=NBMIX+1 IF(CDOOR.EQ.'EXCELL') THEN CALL LCMGET(IPTRK,'STATE-VECTOR',ISTATE) INSB=ISTATE(22) * recover the number of tracks dispached in eack OpenMP core NBATCH=ISTATE(27) IF(NBATCH.EQ.0) NBATCH=1 IF((INSB.GE.1).AND.(ISTATE(7).NE.5)) GO TO 110 ELSE NBATCH=1 ENDIF *---- * COMPUTE THE REDUCED PIJ MATRIX -- NON-VECTORIAL ALGORITHM. *---- ALLOCATE(PGAR(IPIJK*NELPIJ),SGAR(NB1),SGAS(NB1*NANI)) ALLOCATE(ALBP(NALBP,1)) DO 100 IGR=1,NGRP IOFSET=NPSYS(IGR) IF(IOFSET.NE.0) THEN IF(IMPX.GT.10) WRITE(IOUT,'(/25H DOORPV: PROCESSING GROUP,I5, > 6H WITH ,A,1H.)') IGR,CDOOR KPSYS=LCMGIL(JPSYS,IOFSET) IF(LBIHET) CALL LCMSIX(KPSYS,'BIHET',1) CALL LCMGET(KPSYS,'DRAGON-TXSC',SGAR) CALL LCMLEN(KPSYS,'DRAGON-S0XSC',ILONG,ITYLCM) IF(ILONG.GT.NB1*NANI) CALL XABORT('DOORPV: S0XSC OVERFLOW(1).') IF(MOD(ITPIJ,2).EQ.1) THEN CALL LCMGET(KPSYS,'DRAGON-S0XSC',SGAS) ELSE ! avoid scattering reduction SGAS(:NB1*NANI)=0.0 ENDIF IF(NALBP.GT.0) CALL LCMGET(KPSYS,'ALBEDO',ALBP) IF((CDOOR.EQ.'EXCELL').AND.(ISTATE(7).EQ.5)) THEN ! MULTICELL SURFACIC APPROXIMATION IF(ISTATE(10).NE.0) CALL XABORT('DOORPV: TISO EXPECTED.') CALL MUSP(IPTRK,IFTRAK,IMPX,NREG,NBMIX,MAT,VOL,SGAR,SGAS, > NELPIJ,LEAKSW,NBATCH,TITR,NALBP,ALBP,PGAR) ELSE IF(CDOOR.EQ.'EXCELL') THEN IF(IPIJK.EQ.4) THEN NNREG=NREG*NREG NPST=3*NNREG ALLOCATE(PROBKS(NPST)) ELSE NPST=1 ENDIF IF(IFTRAK .NE. 0) REWIND IFTRAK N2PRO= (NREG+NSOUT+1)**2 ALLOCATE(MATALB(-NSOUT:NREG),VOLSUR(-NSOUT:NREG,1), > DPROB(N2PRO,1),DPROBX(N2PRO,1)) CALL EXCELP(IPTRK,IFTRAK,IMPX,NSOUT,NREG,NBMIX,MAT,KNORM,SGAR, > IPIJK,N2PRO,1,NNPSYS(1),NBATCH,TITR,NALBP,ALBP,MATALB,VOLSUR, > DPROB,DPROBX) CALL PIJWIJ(IPTRK,IMPX,NSOUT,NREG,NBMIX,NANI,MAT,VOL,SGAR, > SGAS,NELPIJ,IPIJK,LEAKSW,N2PRO,1,NNPSYS(1),NPST,NALBP,ALBP, > MATALB,VOLSUR,DPROB,DPROBX,PGAR(1),PROBKS) DEALLOCATE(DPROBX,DPROB,VOLSUR,MATALB) IF(IPIJK.EQ.4) THEN CALL LCMPUT(KPSYS,'DRAGON1P*SCT',NNREG,2,PROBKS) CALL LCMPUT(KPSYS,'DRAGON2P*SCT',NNREG,2,PROBKS(NNREG+1)) CALL LCMPUT(KPSYS,'DRAGON3P*SCT',NNREG,2,PROBKS(2*NNREG+1)) DEALLOCATE(PROBKS) ENDIF ELSE IF(CDOOR(:5).EQ.'SYBIL') THEN CALL SYBILP(IPTRK,IMPX,NREG,NBMIX,MAT,VOL,SGAR,SGAS,NELPIJ, > PGAR,LEAKSW) ELSE CALL XABORT('DOORPV: UNKNOWN PIJ DOOR NAMED '//CDOOR//'.') ENDIF *---- * REMOVE LEAKAGE FROM THE SCATTERING-REDUCED CP MATRIX. *---- IF(LNORM) THEN ALLOCATE(PIS(NREG)) CALL XDRNRM(NREG,NBMIX,MAT,VOL,SGAR(2),SGAS(2),PGAR(1),PIS) IF(LEAKSW) CALL LCMPUT(KPSYS,'DRAGON-WIS',NREG,2,PIS) DEALLOCATE(PIS) ENDIF *---- * FORMAT THE REDUCED PIJ MATRIX AS A SQUARE MATRIX. *---- PREG_PTR=LCMARA(NREG*NREG) CALL C_F_POINTER(PREG_PTR,PREG,(/ NREG,NREG /)) DO 15 I=1,NREG FACT=1.0/VOL(I) DO 10 J=1,NREG PREG(I,J)=PGAR(INDPOS(I,J))*FACT 10 CONTINUE 15 CONTINUE CALL LCMPPD(KPSYS,'DRAGON-PCSCT',NREG*NREG,2,PREG_PTR) IF(IPIJK.EQ.4) THEN DO 30 IJKS=1,3 PREG_PTR=LCMARA(NREG*NREG) CALL C_F_POINTER(PREG_PTR,PREG,(/ NREG,NREG /)) DO 25 I=1,NREG FACT=1.0/VOL(I) DO 20 J=1,NREG KS=NELPIJ*IJKS+INDPOS(I,J) PREG(I,J)=PGAR(KS)*FACT 20 CONTINUE 25 CONTINUE CALL LCMPPD(KPSYS,'DRAGON'//CMS(IJKS)//'PCSCT',NREG*NREG,2, > PREG_PTR) 30 CONTINUE ENDIF IF(LBIHET) CALL LCMSIX(KPSYS,' ',2) ENDIF 100 CONTINUE DEALLOCATE(ALBP,SGAS,SGAR,PGAR) GO TO 210 *---- * COMPUTE THE REDUCED PIJ MATRIX -- VECTORIAL ALGORITHM FOR EXCELP. *---- 110 ALLOCATE(SGAR(NB1*NGRP),SGAS(NB1*NANI*NGRP)) ALLOCATE(ALBP(NALBP,NGRP)) DO 120 IGR=1,NGRP IOFSET=NPSYS(IGR) IF(IOFSET.NE.0) THEN KPSYS=LCMGIL(JPSYS,IOFSET) IF(LBIHET) CALL LCMSIX(KPSYS,'BIHET',1) CALL LCMGET(KPSYS,'DRAGON-TXSC',SGAR((IGR-1)*NB1+1)) CALL LCMLEN(KPSYS,'DRAGON-S0XSC',ILONG,ITYLCM) IF(ILONG.GT.NB1*NANI) CALL XABORT('DOORPV: S0XSC OVERFLOW(2).') IF(MOD(ITPIJ,2).EQ.1) THEN CALL LCMGET(KPSYS,'DRAGON-S0XSC',SGAS((IGR-1)*NB1*NANI+1)) ELSE ! avoid scattering reduction SGAS((IGR-1)*NB1*NANI+1:IGR*NB1*NANI)=0.0 ENDIF IF(NALBP.GT.0) CALL LCMGET(KPSYS,'ALBEDO',ALBP(1,IGR)) IF(LBIHET) CALL LCMSIX(KPSYS,' ',2) ENDIF 120 CONTINUE ALLOCATE(PGARG(IPIJK*NELPIJ,NGRP)) IF((ISTATE(7).EQ.4).OR.(INSB.EQ.1)) THEN * --- ALLG VECTORIZATION IF(IPIJK.EQ.4) THEN NNREG=NREG*NREG NPST=3*NNREG ALLOCATE(PROBKS(NPST*NGRP)) ELSE NPST=1 ENDIF IF(IFTRAK .NE. 0) REWIND IFTRAK N2PRO= (NREG+NSOUT+1)**2 ALLOCATE(MATALB(-NSOUT:NREG),VOLSUR(-NSOUT:NREG,NGRP), > DPROB(N2PRO,NGRP),DPROBX(N2PRO,NGRP)) CALL EXCELP(IPTRK,IFTRAK,IMPX,NSOUT,NREG,NBMIX,MAT,KNORM,SGAR, > IPIJK,N2PRO,NGRP,NPSYS(1),NBATCH,TITR,NALBP,ALBP,MATALB,VOLSUR, > DPROB,DPROBX) CALL PIJWIJ(IPTRK,IMPX,NSOUT,NREG,NBMIX,NANI,MAT,VOL,SGAR,SGAS, > NELPIJ,IPIJK,LEAKSW,N2PRO,NGRP,NPSYS(1),NPST,NALBP,ALBP,MATALB, > VOLSUR,DPROB,DPROBX,PGARG(1,1),PROBKS) DEALLOCATE(DPROBX,DPROB,VOLSUR,MATALB) ELSE IF(INSB.EQ.2) THEN * --- XCLL VECTORIZATION IF(IPIJK.NE.1) CALL XABORT('DOORPV: INVALID VALUE OF IPIJK') CALL PIJXL3(IPTRK,IMPX,NGRP,NANI,NBMIX,NPSYS,KNORM,LEAKSW, > SGAR,SGAS,NELPIJ,PGARG) ELSE CALL XABORT('DOORPV: INVALID VALUE OF INSB') ENDIF KPST=0 DO 200 IGR=1,NGRP IOFSET=NPSYS(IGR) IF(IOFSET.NE.0) THEN KPSYS=LCMGIL(JPSYS,IOFSET) IF(LBIHET) CALL LCMSIX(KPSYS,'BIHET',1) IF(IPIJK.EQ.4) THEN CALL LCMPUT(KPSYS,'DRAGON1P*SCT',NNREG,2,PROBKS(KPST+1)) CALL LCMPUT(KPSYS,'DRAGON2P*SCT',NNREG,2,PROBKS(KPST+NNREG+1)) CALL LCMPUT(KPSYS,'DRAGON3P*SCT',NNREG,2, > PROBKS(KPST+2*NNREG+1)) ENDIF *---- * REMOVE LEAKAGE FROM THE SCATTERING-REDUCED CP MATRIX. *---- IF(LNORM) THEN ALLOCATE(PIS(NREG)) CALL XDRNRM(NREG,NBMIX,MAT,VOL,SGAR((IGR-1)*NB1+2), > SGAS((IGR-1)*NB1*NANI+2),PGARG(1,IGR),PIS) IF(LEAKSW) CALL LCMPUT(KPSYS,'DRAGON-WIS',NREG,2,PIS) DEALLOCATE(PIS) ENDIF *---- * FORMAT THE REDUCED PIJ MATRIX AS A SQUARE MATRIX. *---- PREG_PTR=LCMARA(NREG*NREG) CALL C_F_POINTER(PREG_PTR,PREG,(/ NREG,NREG /)) DO 135 I=1,NREG FACT=1.0/VOL(I) DO 130 J=1,NREG PREG(I,J)=PGARG(INDPOS(I,J),IGR)*FACT 130 CONTINUE 135 CONTINUE CALL LCMPPD(KPSYS,'DRAGON-PCSCT',NREG*NREG,2,PREG_PTR) IF(IPIJK.EQ.4) THEN DO 150 IJKS=1,3 PREG_PTR=LCMARA(NREG*NREG) CALL C_F_POINTER(PREG_PTR,PREG,(/ NREG,NREG /)) DO 145 I=1,NREG FACT=1.0/VOL(I) DO 140 J=1,NREG KS=NELPIJ*IJKS+INDPOS(I,J) PREG(I,J)=PGARG(KS,IGR)*FACT 140 CONTINUE 145 CONTINUE CALL LCMPPD(KPSYS,'DRAGON'//CMS(IJKS)//'PCSCT',NREG*NREG,2, > PREG_PTR) 150 CONTINUE ENDIF IF(LBIHET) CALL LCMSIX(KPSYS,' ',2) ENDIF KPST=KPST+NPST 200 CONTINUE IF(IPIJK.EQ.4) DEALLOCATE(PROBKS) DEALLOCATE(ALBP,SGAS,SGAR,PGARG) *---- * DOUBLE HETEROGENEITY TREATMENT *---- 210 IF(LBIHET) THEN NREG=NREGAR NBMIX=NBMIXG DO I=1,NREG MAT(I)=MAT2(I) VOL(I)=VOL2(I) ENDDO DEALLOCATE(MAT2,VOL2) ENDIF RETURN END