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
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
|
*DECK TRIRCA
SUBROUTINE TRIRCA(IPMACR,IPMACP,NGRP,NBMIX,NANI,NW,LDIFF,IL,IPR,
1 RCAT)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Compute the RCAT removal matrix in SPN cases.
*
*Copyright:
* Copyright (C) 2012 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
* IPMACR L_MACROLIB pointer to the unperturbed cross sections.
* IPMACP L_MACROLIB pointer to the perturbed cross sections if
* IPR.gt.0. Equal to IPMACR if IPR=0.
* NGRP number of energy groups.
* NBMIX total number of material mixtures in the macrolib.
* NANI maximum scattering order recovered from tracking and macrolib.
* NW maximum Legendre order (0 or 1) for the total cross sections.
* LDIFF flag set to .true. to use 1/3D as 'NTOT1' cross sections.
* IL scattering Legendre order.
* IPR type of assembly:
* =0: calculation of the system matrices;
* =1: calculation of the derivative of these matrices;
* =2: calculation of the first variation of these matrices;
* =3: identical to IPR=2, but these variation are added to
* unperturbed system matrices.
*
*Parameters: output
* RCAT removal matrix.
*
*-----------------------------------------------------------------------
*
USE GANLIB
*----
* SUBROUTINE ARGUMENTS
*----
TYPE(C_PTR) IPMACR,IPMACP
INTEGER NGRP,NBMIX,NANI,NW,IL,IPR
LOGICAL LDIFF
DOUBLE PRECISION RCAT(NGRP,NGRP,NBMIX)
*----
* LOCAL VARIABLES
*----
TYPE(C_PTR) JPMACR,JPMACP,KPMACR,KPMACP
CHARACTER TEXT12*12,CM*2,HSMG*131
DOUBLE PRECISION OTH
INTEGER, ALLOCATABLE, DIMENSION(:) :: IJJ,NJJ,IPOS
REAL, ALLOCATABLE, DIMENSION(:) :: WORK
REAL, ALLOCATABLE, DIMENSION(:,:) :: SGD
PARAMETER(OTH=1.0D0/3.0D0)
*----
* SCRATCH STORAGE ALLOCATION
*----
ALLOCATE(IJJ(NBMIX),NJJ(NBMIX),IPOS(NBMIX))
ALLOCATE(SGD(NBMIX,3),WORK(NBMIX*NGRP))
*
JPMACR=LCMGID(IPMACR,'GROUP')
JPMACP=LCMGID(IPMACP,'GROUP')
WRITE(CM,'(I2.2)') IL-1
RCAT(:NGRP,:NGRP,:NBMIX)=0.0D0
DO 100 IGR=1,NGRP
* PROCESS SECONDARY GROUP IGR.
KPMACP=LCMGIL(JPMACP,IGR)
SGD(:NBMIX,1)=0.0
CALL LCMLEN(KPMACP,'SIGW'//CM,LENGT,ITYLCM)
IF((LENGT.GT.0).AND.(IL.LE.NANI)) THEN
IF(LENGT.GT.NBMIX) CALL XABORT('TRIRCA: INVALID LENGTH FOR'
1 //' SIGW'//CM//' CROSS SECTIONS.')
CALL LCMGET(KPMACP,'SIGW'//CM,SGD(1,1))
ENDIF
WRITE(TEXT12,'(4HNTOT,I1)') MIN(IL-1,9)
CALL LCMLEN(KPMACP,TEXT12,LENGT,ITYLCM)
CALL LCMLEN(KPMACP,'NTOT1',LENGT1,ITYLCM)
IF((IL.EQ.1).AND.(LENGT.NE.NBMIX)) CALL XABORT('TRIRCA: NO NTOT0'
1 //' CROSS SECTIONS.')
IF(MOD(IL-1,2).EQ.0) THEN
* macroscopic total cross section in even-parity equations.
IF(LENGT.EQ.NBMIX) THEN
CALL LCMGET(KPMACP,TEXT12,SGD(1,2))
ELSE
CALL LCMGET(KPMACP,'NTOT0',SGD(1,2))
ENDIF
DO 10 IBM=1,NBMIX
IF((SGD(IBM,2)-SGD(IBM,1).LT.0.0).AND.(IPR.EQ.0)) THEN
WRITE(HSMG,'(28HTRIRCA: NEGATIVE XS IN GROUP,I5)') IGR
CALL XABORT(HSMG)
ENDIF
RCAT(IGR,IGR,IBM)=SGD(IBM,2)-SGD(IBM,1)
10 CONTINUE
ELSE
* macroscopic total cross section in odd-parity equations.
IF(LDIFF) THEN
CALL LCMLEN(KPMACP,'DIFF',LENGT,ITYLCM)
IF(LENGT.EQ.0) CALL XABORT('TRIRCA: DIFFUSION COEFFICIENTS'
1 //' EXPECTED IN THE MACROLIB.')
IF(LENGT.GT.NBMIX) CALL XABORT('TRIRCA: INVALID LENGTH FOR'
1 //' DIFFUSION COEFFICIENTS.')
CALL LCMGET(KPMACP,'DIFF',SGD(1,2))
IF(IPR.EQ.0) THEN
DO 20 IBM=1,NBMIX
RCAT(IGR,IGR,IBM)=OTH/SGD(IBM,2)
20 CONTINUE
ELSE IF(IPR.EQ.1) THEN
KPMACR=LCMGIL(JPMACR,IGR)
CALL LCMGET(KPMACR,'DIFF',SGD(1,3))
DO 30 IBM=1,NBMIX
RCAT(IGR,IGR,IBM)=-OTH*SGD(IBM,2)/SGD(IBM,3)**2
30 CONTINUE
ELSE IF(IPR.EQ.2) THEN
KPMACR=LCMGIL(JPMACR,IGR)
CALL LCMGET(KPMACR,'DIFF',SGD(1,3))
DO 40 IBM=1,NBMIX
RCAT(IGR,IGR,IBM)=OTH/(SGD(IBM,2)+SGD(IBM,3))
1 -OTH/SGD(IBM,3)
40 CONTINUE
ELSE IF(IPR.EQ.3) THEN
KPMACR=LCMGIL(JPMACR,IGR)
CALL LCMGET(KPMACR,'DIFF',SGD(1,3))
DO 50 IBM=1,NBMIX
RCAT(IGR,IGR,IBM)=OTH/(SGD(IBM,2)+SGD(IBM,3))
50 CONTINUE
ENDIF
GO TO 100
ELSE
IF(NW.EQ.0) THEN
CALL LCMGET(KPMACP,'NTOT0',SGD(1,2))
ELSE IF(LENGT.EQ.NBMIX) THEN
CALL LCMGET(KPMACP,TEXT12,SGD(1,2))
ELSE IF((NW.GE.1).AND.(LENGT1.EQ.NBMIX)) THEN
CALL LCMGET(KPMACP,'NTOT1',SGD(1,2))
ELSE
CALL LCMGET(KPMACP,'NTOT0',SGD(1,2))
ENDIF
DO 60 IBM=1,NBMIX
RCAT(IGR,IGR,IBM)=SGD(IBM,2)-SGD(IBM,1)
60 CONTINUE
ENDIF
IF(IPR.EQ.0) THEN
DO 65 IBM=1,NBMIX
IF(RCAT(IGR,IGR,IBM).LT.0.0) THEN
WRITE(HSMG,'(39HTRIRCA: INVALID CROSS-SECTION DATA (IL=,
1 I3,2H).)') IL
CALL XABORT(HSMG)
ENDIF
65 CONTINUE
ENDIF
ENDIF
CALL LCMLEN(KPMACP,'NJJS'//CM,LENGT,ITYLCM)
IF(LENGT.GT.NBMIX) CALL XABORT('TRIRCA: INVALID LENGTH FOR NJJS'
1 //CM//' INFORMATION.')
IF((LENGT.GT.0).AND.(IL.LE.NANI)) THEN
CALL LCMGET(KPMACP,'NJJS'//CM,NJJ)
CALL LCMGET(KPMACP,'IJJS'//CM,IJJ)
IGMIN=IGR
IGMAX=IGR
DO 70 IBM=1,NBMIX
IGMIN=MIN(IGMIN,IJJ(IBM)-NJJ(IBM)+1)
IGMAX=MAX(IGMAX,IJJ(IBM))
70 CONTINUE
CALL LCMGET(KPMACP,'IPOS'//CM,IPOS)
CALL LCMGET(KPMACP,'SCAT'//CM,WORK)
DO 90 JGR=IGMAX,IGMIN,-1
IF(JGR.EQ.IGR) GO TO 90
DO 80 IBM=1,NBMIX
IF((JGR.GT.IJJ(IBM)-NJJ(IBM)).AND.(JGR.LE.IJJ(IBM))) THEN
RCAT(IGR,JGR,IBM)=-WORK(IPOS(IBM)+IJJ(IBM)-JGR)
ENDIF
80 CONTINUE
90 CONTINUE
ENDIF
100 CONTINUE
*----
* SCRATCH STORAGE DEALLOCATION
*----
DEALLOCATE(WORK,SGD)
DEALLOCATE(IPOS,NJJ,IJJ)
RETURN
END
|