1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
|
*DECK PIJRGL
SUBROUTINE PIJRGL(IPRT,NREG,NSOUT,SIGTAL,PROB)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Gelbard normalization of collision probs (CP).
*
*Copyright:
* Copyright (C) 1991 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, G. Marleau
*
*Parameters: input
* IPRT print level.
* NREG number of zones for geometry.
* NSOUT number of surfaces for geometry.
* SIGTAL albedo-sigt vector.
*
*Parameters: input/output
* PROB CP matrix for all types.
*
*References:
* R. Roy and G. Marleau,
* Normalization Techniques for CP Matrices,
* CONF/PHYSOR-90, Marseille/France, V 2, P IX-40 (1990).
*-----------------------------------------------------------------------
*
IMPLICIT NONE
INTEGER IPRT,NREG,NSOUT,IUNOUT,IPRINT,
> IPRB,IPRF,IUNK,JUNK,IVOL,JVOL,IR,JR,
> NSURM,NSURC,NVOLM,NVOLC,IP
PARAMETER (IUNOUT=6, IPRINT=4)
REAL SIGTAL(-NSOUT:NREG)
DOUBLE PRECISION PROB(*),RBARRE,GBARRE
DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: RI
C----
C SCRATCH STORAGE ALLOCATION
C----
ALLOCATE(RI(-NSOUT:NREG))
C
RBARRE=0.0
GBARRE=0.0
IPRB= 0
IUNK= 0
IVOL= NSOUT*(NSOUT+1)/2
C
C COMPUTE R-SUB(I) FACTORS AND: GBARRE, RBARRE
DO 30 IR=-NSOUT, NREG
IUNK= IUNK+1
RI(IR)=0.0
DO 10 JR=-NSOUT, IR
IPRB= IPRB+1
IF( JR.LT.0.OR.SIGTAL(JR).GT.0.0 )THEN
RI(IR)=RI(IR)+PROB(IPRB)
ENDIF
10 CONTINUE
IPRF= IPRB
JUNK= IUNK
DO 20 JR= IR+1 , NREG
IPRF= IPRF+JUNK
JUNK= JUNK+1
IF( JR.LT.0.OR.SIGTAL(JR).GT.0.0 )THEN
RI(IR)=RI(IR)+PROB(IPRF)
ENDIF
20 CONTINUE
IF( IR.LT.0 )THEN
IVOL= IVOL+1
RI(IR)= PROB(IVOL)-RI(IR)
GBARRE= GBARRE+PROB(IVOL)
RBARRE= RBARRE+RI(IR)
ELSEIF( IR.GT.0 )THEN
IVOL= IVOL+IUNK-1
RI(IR)= PROB(IVOL)-RI(IR)
IF( SIGTAL(IR).GT.0.0 )THEN
GBARRE= GBARRE+PROB(IVOL)
RBARRE= RBARRE+RI(IR)
ENDIF
ELSE
IVOL= IVOL+1
RI(IR)=0.0
ENDIF
30 CONTINUE
GBARRE=1.0/GBARRE
RBARRE=RBARRE*GBARRE
C
C RENORMALIZE PROB MATRIX
IVOL= NSOUT*(NSOUT+1)/2
IPRB= 0
IUNK= 0
DO 210 IR = -NSOUT, NREG
IF( IR.LE.0 )THEN
IVOL= IVOL+1
ELSE
IVOL= IVOL+IUNK
ENDIF
IUNK= IUNK+1
JVOL= NSOUT*(NSOUT+1)/2
JUNK= 0
DO 200 JR= -NSOUT, IR
IF( JR.LE.0 )THEN
JVOL= JVOL+1
ELSE
JVOL= JVOL+JUNK
ENDIF
JUNK= JUNK+1
IPRB= IPRB+1
IF( IR.NE.0.AND.JR.NE.0 )THEN
PROB(IPRB)= PROB(IPRB)+(PROB(JVOL)*RI(IR)
> +PROB(IVOL)*RI(JR)-PROB(IVOL)*PROB(JVOL)*RBARRE)*GBARRE
ENDIF
200 CONTINUE
210 CONTINUE
C
C PRINT IF REQUESTED
IF( IPRT.GE.IPRINT )THEN
WRITE(IUNOUT,'(19H0 GLOBAL FACTORS: ,
> 8H RBARRE=,1P,F11.5,5X,7HGBARRE=,F11.5)')
> RBARRE, GBARRE
WRITE(IUNOUT,'(30H0 SURFACE ADJUSTMENT FACTORS /)')
NSURC = -1
DO 300 IP = 1, (9 +NSOUT) / 10
NSURM= MAX( -NSOUT, NSURC-9 )
WRITE(IUNOUT,'(10X,10( A5, I6)/)')
> (' SUR ',-IR,IR= NSURC, NSURM, -1)
WRITE(IUNOUT,'(10H R-SUB(I) ,10F11.5)')
> (RI(IR),IR=NSURC,NSURM,-1)
NSURC = NSURC - 10
300 CONTINUE
WRITE(IUNOUT,'(30H0 VOLUME ADJUSTMENT FACTORS /)')
NVOLC = 1
DO 310 IP = 1, (9 + NREG) / 10
NVOLM= MIN( NREG, NVOLC+9 )
WRITE(IUNOUT,'(10X,10( A5 , I6)/)')
> (' VOL ',IR,IR=NVOLC,NVOLM, 1)
WRITE(IUNOUT,'(10H R-SUB(I) ,10F11.5)')
> (RI(IR),IR=NVOLC,NVOLM, 1)
NVOLC = NVOLC + 10
310 CONTINUE
ENDIF
C----
C SCRATCH STORAGE DEALLOCATION
C----
DEALLOCATE(RI)
C
RETURN
END
|