diff options
Diffstat (limited to 'Dragon/src/PIJCMP.f')
| -rw-r--r-- | Dragon/src/PIJCMP.f | 121 |
1 files changed, 121 insertions, 0 deletions
diff --git a/Dragon/src/PIJCMP.f b/Dragon/src/PIJCMP.f new file mode 100644 index 0000000..f62fb2e --- /dev/null +++ b/Dragon/src/PIJCMP.f @@ -0,0 +1,121 @@ +*DECK PIJCMP + SUBROUTINE PIJCMP(NREG,NSOUT,NCOR,DPR,VOLSUR,LPIJK,PROB) +* +*----------------------------------------------------------------------- +* +*Purpose: +* Compression of PIJ matrices in symmetric format. +* +*Copyright: +* Copyright (C) 1994 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): R. Roy +* +*Parameters: input +* NREG number of zones for geometry. +* NSOUT number of surfaces for geometry. +* NCOR maximum number of corners. +* DPR collision probabilities. +* VOLSUR volumes. +* LPIJK pijk flag. +* +*Parameters: output +* PROB compress probability matrix. +* +*Comments: +* Format of compress probability matrix +* NPLEN=(NREG+NSOUT+2)*(NREG+NSOUT+1)/2 +* IND(I,J)=MAX(I+NSOUT+1,J+NSOUT+1) +* *(MAX(I+NSOUT+1,J+NSOUT+1)-1)/2 +* +MIN(I+NSOUT+1,J+NSOUT+1) +* IS=-NSOUT,-1; JS=-NSOUT,IS; I=IND(IS,JS) +* PROB(I)=VOLSUR(IS)*PSS(IS,JS) +* IV=1,NREG; JS=-NSOUT,-1; I=IND(IV,JS) +* SIGT(IV).GT.0.0 +* PROB(I)=SIGT(IV)*VOLSUR(IV)*PVS(IV,JS) +* SIGT(IV).EQ.0.0 +* PROB(I)=VOLSUR(IV)*PVS(IV,JS) +* IV=1,NREG; JV=1,IV; I=IND(IV,JV) +* SIGT(IV).GT.0.0 AND SIGT(JV).GT.0.0 +* PROB(I)=SIGT(IV)*SIGT(JV)*VOLSUR(IV)*PVV(IV,JV) +* SIGT(IV).GT.0.0 AND SIGT(JV).EQ.0.0 +* PROB(I)=SIGT(IV)*VOLSUR(IV)*PVV(IV,JV) +* SIGT(IV).EQ.0.0 AND SIGT(JV).GT.0.0 +* PROB(I)=SIGT(JV)*VOLSUR(IV)*PVV(IV,JV) +* SIGT(IV).EQ.0.0 AND SIGT(JV).EQ.0.0 +* PROB(I)=VOLSUR(IV)*PVV(IV,JV) +* +*----------------------------------------------------------------------- +* + IMPLICIT NONE +*---- +* VARIABLES +*---- + INTEGER NREG,NSOUT,NCOR,IPR,IL,JL,IVOL,IUN + DOUBLE PRECISION DPR(-NSOUT:NREG,-NSOUT:NREG),PROB(*),ZCOR,ZCOR1 + REAL VOLSUR(-NSOUT:NREG),COEF + LOGICAL LPIJK + DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: WORK +*---- +* SCRATCH STORAGE ALLOCATION +*---- + ALLOCATE(WORK((NREG+NSOUT+2)*(NREG+NSOUT+1)/2)) +*---- +* SYMMETRIZE AND STORE IN PROB +*---- + IPR= 0 + DO 150 IL = -NSOUT, NREG + DO 160 JL = -NSOUT, IL + IPR= IPR+1 + WORK(IPR)= DPR(IL,JL) + DPR(JL,IL) + 160 CONTINUE + 150 CONTINUE + IF( NCOR.EQ.1 )THEN + IPR= 0 + DO 250 IL = -NSOUT, NREG + DO 260 JL = -NSOUT, IL + IPR= IPR+1 + PROB(IPR)= WORK(IPR) + 260 CONTINUE + 250 CONTINUE + ELSE + IPR= 0 + ZCOR1= 1.0D0/DBLE(NCOR) + ZCOR= 1.0D0/DBLE(NCOR*NCOR) + DO 251 IL = -NSOUT, NREG + IF( IL.GT.0 ) ZCOR= ZCOR1 + DO 261 JL = -NSOUT, IL + IPR= IPR+1 + IF( JL.GT.0 ) ZCOR= 1.0D0 + PROB(IPR)= ZCOR * WORK(IPR) + 261 CONTINUE + 251 CONTINUE + ENDIF +*---- +* CHARGE VOLUMES IN THE PROB MATRIX +*---- + COEF=1.0 + IVOL= NSOUT*(NSOUT+1)/2 + DO 300 IUN= -NSOUT, NREG + IF( IUN.LE.0 )THEN + IVOL= IVOL+1 + IF(LPIJK) COEF= 3./4. + ELSE + IVOL= IVOL+NSOUT+IUN + IF(LPIJK) COEF= 2./3. + ENDIF + IF( PROB(IVOL).NE.0.0 )THEN + CALL XABORT( 'PIJCMP: UNEXPECTED VALUE IN PROB MATRIX' ) + ENDIF + PROB(IVOL) = VOLSUR(IUN)*COEF + 300 CONTINUE +*---- +* SCRATCH STORAGE DEALLOCATION +*---- + DEALLOCATE(WORK) + RETURN + END |