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
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
|
*DECK AUTFLU
SUBROUTINE AUTFLU(IPTRK,IPLIB,IPLI0,IFTRAK,NREG,NUN,NBMIX,NBISO,
1 NIRES,MAT,VOL,KEYFLX,CDOOR,LEAKSW,IMPX,DEN,MIX,IAPT,IPHASE,NGRP,
2 IGRMIN,IGRRES,IGRMAX,DIL,TITR,IALTER,DELI,LBIN,NBIN,EBIN,MAXTRA,
3 ISEED,ITRANC,UUU,FUNKNO,SIGT,SIGS,SIGS1,SIGF,SIGGAR,MASKG)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Compute the flux in the Autolib fine groups using the Autosecol
* method.
*
*Copyright:
* Copyright (C) 2023 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
* IPTRK pointer to the tracking (L_TRACK signature).
* IPLIB pointer to the internal microscopic cross section library
* with subgroups (L_LIBRARY signature).
* IPLI0 pointer to the internal microscopic cross section library
* builded by the self-shielding module.
* IFTRAK file unit number used to store the tracks.
* NREG number of regions.
* NUN number of unknowns per energy group and band.
* NBMIX number of mixtures in the internal library.
* NBISO number of isotopes.
* NIRES number of correlated resonant isotopes.
* MAT index-number of the mixture type assigned to each volume.
* VOL volumes.
* KEYFLX pointers of fluxes in unknown vector.
* CDOOR name of the geometry/solution operator.
* LEAKSW leakage flag (LEAKSW=.true. if neutron leakage through
* external boundary is present).
* IMPX print flag (equal to zero for no print).
* DEN density of each isotope.
* MIX mix number of each isotope (can be zero).
* IAPT resonant isotope index associated with isotope I. Mixed
* moderator if IAPT(I)=NIRES+1. Out-of-fuel isotope if
* IAPT(I)=0.
* IPHASE type of flux solution (=1 use a native flux solution door;
* =2 use collision probabilities).
* NGRP number of energy groups.
* IGRMIN first group where the self-shielding is applied.
* IGRRES first resolved energy group.
* IGRMAX most thermal group where the self-shielding is applied.
* DIL microscopic dilution cross section of each isotope.
* TITR title.
* IALTER type of elastic slowing-down kernel (=0: use exact kernel;
* =1: use an approximate kernel for the resonant isotopes).
* DELI elementary lethargy width used by the elastic kernel.
* LBIN total number of fine energy groups in the Autolib.
* NBIN number of fine energy groups in each coarse energy group.
* EBIN energy limits of the Autolib fine groups.
* MAXTRA maximum number of down-scattering terms.
* ISEED the seed for the generation of random numbers in the
* unresolved energy domain.
* ITRANC type of transport correction.
*
*Parameters: output
* UUU lethargy limits of the Autolib fine groups.
* FUNKNO flux in the Autolib fine groups.
* SIGT total microscopic x-s.
* SIGS P0 scattering microscopic x-s.
* SIGS1 P1 scattering microscopic x-s.
* SIGF nu*fission microscopic x-s.
* SIGGAR macroscopic x-s of the non-resonant isotopes in each mixture:
* (*,*,*,1) total; (*,*,*,2) transport correction;
* (*,*,*,3) P0 scattering.
* MASKG energy group mask pointing on self-shielded groups.
*
*-----------------------------------------------------------------------
*
USE GANLIB
*----
* SUBROUTINE ARGUMENTS
*----
TYPE(C_PTR) IPTRK,IPLIB,IPLI0
INTEGER IFTRAK,NREG,NUN,NBMIX,NBISO,NIRES,MAT(NREG),KEYFLX(NREG),
1 IMPX,MIX(NBISO),IAPT(NBISO),IPHASE,NGRP,IGRMIN,IGRRES,IGRMAX,
2 IALTER,LBIN,NBIN(NGRP),MAXTRA,ISEED,ITRANC
REAL VOL(NREG),DEN(NBISO),DIL(NBISO),DELI,EBIN(LBIN+1),
1 UUU(LBIN+1),FUNKNO(NUN,LBIN),SIGT(LBIN,NBISO),SIGS(LBIN,NBISO),
2 SIGS1(LBIN,NBISO),SIGF(LBIN,NBISO),SIGGAR(NBMIX,0:NIRES,NGRP,3)
LOGICAL LEAKSW,MASKG(NGRP)
CHARACTER CDOOR*12,TITR*72
*----
* LOCAL VARIABLES
*----
TYPE(C_PTR) IPP,KPLIB,KPLI0,IPSYS
LOGICAL LLIB
DOUBLE PRECISION DUUU
CHARACTER HNAMIS*12,HNABIS*12,HSMG*131
*----
* ALLOCATABLE ARRAYS
*----
TYPE(C_PTR), ALLOCATABLE, DIMENSION(:) :: IPISO1,IPISO2
INTEGER, ALLOCATABLE, DIMENSION(:) :: NJJ,IJJ,NEXT,III
REAL, ALLOCATABLE, DIMENSION(:) :: GA1,DELTAU,SGAR,PRI
REAL, ALLOCATABLE, DIMENSION(:,:) :: SIGINF,SCAT,GA2,CONC
LOGICAL, ALLOCATABLE, DIMENSION(:) :: MASKH
*----
* SCRATCH STORAGE ALLOCATION
*----
ALLOCATE(GA1(NGRP),GA2(NGRP,NGRP),CONC(NBMIX,NBISO),DELTAU(NGRP),
1 MASKH(NGRP),SIGINF(NGRP,3))
ALLOCATE(IPISO1(NBISO),IPISO2(NBISO))
ALLOCATE(NJJ(NGRP),IJJ(NGRP),SGAR(NGRP**2),SCAT(NGRP,NGRP))
ALLOCATE(PRI(MAXTRA),NEXT(NBISO),III(NBISO+1))
NEXT(:NBISO)=0
*
CALL KDRCPU(TK1)
SIGGAR(:NBMIX,0:NIRES,:NGRP,:3)=0.0
CALL LIBIPS(IPLIB,NBISO,IPISO1)
CALL LIBIPS(IPLI0,NBISO,IPISO2)
CALL LCMGET(IPLI0,'DELTAU',DELTAU)
III(1)=1
DO 120 ISO=1,NBISO
SIGT(:LBIN,ISO)=0.0
SIGS(:LBIN,ISO)=0.0
SIGS1(:LBIN,ISO)=0.0
SIGF(:LBIN,ISO)=0.0
IF(IMPX.GT.1) WRITE(6,'(/32H AUTFLU: RECOVER XS FOR ISOTOPE=,I5)')
1 ISO
IBM=MIX(ISO)
DO 10 NRE=1,NREG
IF(MAT(NRE).EQ.IBM) GO TO 20
10 CONTINUE
III(ISO+1)=III(ISO)
GO TO 120
20 KPLIB=IPISO1(ISO) ! infinite dilution isotope
KPLI0=IPISO2(ISO) ! self-shielded isotope
CALL LCMGTC(KPLIB,'ALIAS',12,HNAMIS)
IRES=IAPT(ISO)
JRES=IRES
IF(IRES.EQ.NIRES+1) JRES=0
DENN=DEN(ISO)
CALL LCMGET(KPLIB,'AWR',AWR)
CALL LCMLEN(KPLIB,'BIN-NFS',LENGT,ITYLCM)
IF((IRES.GT.0).AND.(IRES.LE.NIRES).AND.(LENGT.GT.0)) THEN
! resonant isotope
IF(IMPX.GT.2) WRITE(6,'(26H AUTFLU: PROCESS AUTOLIB '',A12,
1 1H'')') HNAMIS
IF(LBIN.EQ.0) CALL XABORT('AUTFLU: MISSING AUTOLIB DATA.')
DUUU=0.0D0
UUU(1)=0.0
DO 40 IGR=1,LBIN
DUUU=DUUU+LOG(EBIN(IGR)/EBIN(IGR+1))
UUU(IGR+1)=REAL(DUUU)
40 CONTINUE
! recover unresolved xs values
LLL=0
IF(IGRRES.GT.IGRMIN) THEN
CALL LCMLEN(KPLIB,'NUSIGF',N10,ITYLCM)
CALL LCMLEN(KPLIB,'SIGS00',N12,ITYLCM)
SIGINF(:NGRP,:3)=0.0
CALL LCMGET(KPLIB,'NTOT0',SIGINF(1,1))
IF(N10.GT.0) CALL LCMGET(KPLIB,'NUSIGF',SIGINF(1,2))
IF(N12.GT.0) CALL LCMGET(KPLIB,'SIGS00',SIGINF(1,3))
CALL AUTTAB(KPLIB,HNAMIS,IGRMIN,IGRRES,NGRP,LBIN,NBIN,UUU,
1 ISEED,SIGINF,LLL,SIGT(1,ISO),SIGS(1,ISO),SIGF(1,ISO))
ENDIF
! recover resolved xs values
IF(LLL.LT.LBIN) THEN
CALL LCMLEN(KPLIB,'BIN-NTOT0',LENG,ITYLCM)
IF(LENG.GT.LBIN) CALL XABORT('AUTFLU: LBIN OVERFLOW.')
CALL LCMGET(KPLIB,'BIN-NTOT0',SIGT(LLL+1,ISO))
CALL LCMGET(KPLIB,'BIN-SIGS00',SIGS(LLL+1,ISO))
CALL LCMLEN(KPLIB,'BIN-NUSIGF',ILEN,ITYLCM)
IF(ILEN.GT.0) CALL LCMGET(KPLIB,'BIN-NUSIGF',SIGF(LLL+1,ISO))
ENDIF
*----
* ELASTIC SCATTERING INFORMATION.
*----
MAXIII=MAXTRA-III(ISO)+1
CALL LIBPRI(MAXIII,DELI,AWR,IALTER,0,NEXT(ISO),PRI(III(ISO)))
ELSE
! infinite dilution isotope
IF(C_ASSOCIATED(KPLI0)) THEN
CALL LCMLEN(KPLI0,'NTOT0',ILEN0,ITYLCM)
IF(ILEN0.NE.0) THEN
LLIB=.FALSE.
IPP=KPLI0
ELSE
LLIB=.TRUE.
IPP=IPISO1(ISO) ! set ISO-th isotope
ENDIF
ELSE
LLIB=.TRUE.
IPP=IPISO1(ISO) ! set ISO-th isotope
ENDIF
IF(LLIB.AND.(.NOT.C_ASSOCIATED(IPP))) THEN
WRITE(HSMG,'(18H AUTFLU: ISOTOPE '',A12,7H'' (ISO=,I8,5H) IS ,
1 39HNOT AVAILABLE IN THE ORIGINAL MICROLIB.)') HNAMIS,ISO
CALL XABORT(HSMG)
ENDIF
IF((.NOT.LLIB).AND.(IMPX.GT.2)) THEN
CALL LCMGTC(KPLI0,'ALIAS',12,HNABIS)
WRITE(6,'(26H AUTFLU: RECOVER ISOTOPE '',A12,11H'' FROM THE ,
1 12HNEW LIBRARY.)') HNABIS
ELSE IF(LLIB.AND.(IMPX.GT.2)) THEN
WRITE(6,'(26H AUTFLU: RECOVER ISOTOPE '',A12,11H'' FROM THE ,
1 12HOLD LIBRARY.)') HNAMIS
ENDIF
NSCAR=0
SCAT(:NGRP,:NGRP)=0.0
CALL LCMLEN(IPP,'NTOT0',N13,ITYLCM)
IF(N13.EQ.0) THEN
CALL LCMLIB(IPP)
CALL XABORT('AUTFLU: NO INFINITE-DILUTION TOTAL XS.')
ELSE IF(N13.GT.LBIN) THEN
CALL XABORT('AUTFLU: LBIN OVERFLOW.')
ELSE IF(N13.NE.NGRP) THEN
CALL XABORT('AUTFLU: INVALID X-SECTIONS.')
ENDIF
CALL LCMGET(IPP,'NTOT0',SIGT(1,ISO))
CALL LCMLEN(IPP,'NUSIGF',N10,ITYLCM)
IF(N10.GT.0) CALL LCMGET(IPP,'NUSIGF',SIGF(1,ISO))
CALL LCMLEN(IPP,'SIGS00',N12,ITYLCM)
IF(N12.GT.0) THEN
CALL LCMGET(IPP,'SIGS00',SIGS(1,ISO))
CALL LCMLEN(IPP,'SIGS01',N14,ITYLCM)
IF(N14.GT.0) THEN
CALL LCMGET(IPP,'SIGS01',SIGS1(1,ISO))
ELSE
SIGS1(:NGRP,ISO)=0.0
ENDIF
ELSE
CALL LCMGET(IPP,'SCAT00',SGAR)
CALL LCMGET(IPP,'NJJS00',NJJ)
CALL LCMGET(IPP,'IJJS00',IJJ)
IGAR1=0
DO IG2=1,NGRP
DO IG1=IJJ(IG2),IJJ(IG2)-NJJ(IG2)+1,-1
IGAR1=IGAR1+1
SCAT(IG1,IG2)=SGAR(IGAR1)
ENDDO
ENDDO
DO IG1=1,NGRP
SUMSC=0.0D0
DO IG2=1,NGRP
SUMSC=SUMSC+SCAT(IG1,IG2)
ENDDO
SIGS(IG1,ISO)=REAL(SUMSC)
ENDDO
CALL LCMLEN(IPP,'SCAT01',N87,ITYLCM)
IF(N87.GT.0) THEN
CALL LCMGET(IPP,'SCAT01',SGAR)
CALL LCMGET(IPP,'NJJS01',NJJ)
CALL LCMGET(IPP,'IJJS01',IJJ)
IGAR1=0
DO IG2=1,NGRP
DO IG1=IJJ(IG2),IJJ(IG2)-NJJ(IG2)+1,-1
IGAR1=IGAR1+1
SCAT(IG1,IG2)=SGAR(IGAR1)
ENDDO
ENDDO
DO IG1=1,NGRP
SUMSC=0.0D0
DO IG2=1,NGRP
SUMSC=SUMSC+SCAT(IG1,IG2)
ENDDO
SIGS1(IG1,ISO)=REAL(SUMSC)
ENDDO
ENDIF
ENDIF
! compute SIGGAR used by SPH equivalence
IF((DENN.NE.0.0).AND.(IBM.NE.0).AND.(JRES.EQ.0)) THEN
DO 70 IG1=1,NGRP
SIGGAR(IBM,JRES,IG1,1)=SIGGAR(IBM,JRES,IG1,1)+DENN*
1 SIGT(IG1,ISO)
SIGGAR(IBM,JRES,IG1,3)=SIGGAR(IBM,JRES,IG1,3)+DENN*
1 SIGS(IG1,ISO)
70 CONTINUE
CALL LCMLEN(IPP,'TRANC',LENGT,ITYLCM)
IF(LENGT.GT.0) THEN
CALL LCMGET(IPP,'TRANC',GA1)
DO 80 IG1=1,NGRP
SIGGAR(IBM,JRES,IG1,2)=SIGGAR(IBM,JRES,IG1,2)+DENN*GA1(IG1)
80 CONTINUE
ENDIF
ENDIF
! expand xs of non-resonant isotopes
CALL AUTPRD(NGRP,LBIN,NBIN,SIGS(1,ISO))
CALL AUTPRD(NGRP,LBIN,NBIN,SIGT(1,ISO))
CALL AUTPRD(NGRP,LBIN,NBIN,SIGF(1,ISO))
CALL AUTPRD(NGRP,LBIN,NBIN,SIGS1(1,ISO))
ENDIF
III(ISO+1)=III(ISO)+NEXT(ISO)
*----
* CROSS SECTION EDITION.
*----
IF(IMPX.GT.7) THEN
CALL LCMGTC(KPLIB,'ALIAS',12,HNAMIS)
WRITE(6,540) HNAMIS,(K,SIGS(K,ISO),SIGT(K,ISO),SIGF(K,ISO),
1 SIGS1(K,ISO),K=1,LBIN)
I2=III(ISO+1)-1
WRITE(6,550) III(ISO),I2,(PRI(K),K=III(ISO),I2)
ENDIF
120 CONTINUE
CALL KDRCPU(TK2)
IF(IMPX.GT.1) WRITE(6,'(/36H AUTFLU: CPU TIME SPENT TO RECOVER A,
1 33HUTOLIB AND INFINITE-DILUTION XS =,F8.1,8H SECOND./)') TK2-TK1
*
CALL KDRCPU(TK1)
TK4=0.0
TK5=0.0
ICPIJ=0
*----
* SET THE NUMBER DENSITIES.
*----
CONC(:NBMIX,:NBISO)=0.0
DO 130 ISO=1,NBISO
IBM=MIX(ISO)
IF(IBM.LE.0) CYCLE
CONC(IBM,ISO)=DEN(ISO)
130 CONTINUE
*----
* DETERMINE WHICH GROUPS ARE SELF-SHIELDED.
*----
MASKG(:NGRP)=.FALSE.
MASKG(IGRMIN:IGRMAX)=.TRUE.
*----
* COMPUTE THE NEUTRON FLUX.
*----
CALL KDRCPU(TKA)
CALL KDRCPU(TKA)
DO 330 IG1=1,NGRP
ICPIJ=ICPIJ+NBIN(IG1)
330 CONTINUE
CALL LCMOP(IPSYS,'**SYSTEM**',0,1,0)
KNORM=1
IMPY=MAX(0,IMPX-3)
IF(IPHASE.EQ.1) THEN
* USE A NATIVE DOOR.
CALL AUTIT2(IPTRK,IFTRAK,IPSYS,MAXTRA,KNORM,NUN,LBIN,NREG,
1 NBMIX,NBISO,MAT,VOL,KEYFLX,NIRES,IAPT,CDOOR,LEAKSW,TITR,IMPY,
2 CONC,SIGS,SIGT,SIGS1,DIL,PRI,UUU,DELI,ITRANC,NEXT,III,FUNKNO)
ELSE IF(IPHASE.EQ.2) THEN
* USE A COLLISION PROBABILITY DOOR.
CALL AUTIT1(IPTRK,IFTRAK,IPSYS,MAXTRA,KNORM,LBIN,NREG,
1 NBMIX,NBISO,MAT,VOL,NIRES,IAPT,CDOOR,LEAKSW,TITR,IMPY,CONC,
2 SIGS,SIGT,SIGS1,DIL,PRI,UUU,DELI,ITRANC,NEXT,III,FUNKNO)
ENDIF
CALL LCMCL(IPSYS,2)
CALL KDRCPU(TKB)
TK4=TK4+(TKB-TKA)
* ***************************************************************
CALL LCMVAL(IPLI0,' ')
*----
* RESET MASKG FOR SPH CALCULATION IN SMALL LETHARGY WIDTH GROUPS.
*----
DO 380 IG1=1,NGRP
IF(MASKG(IG1)) THEN
IF(DELTAU(IG1).GT.0.1) MASKG(IG1)=.TRUE.
ENDIF
380 CONTINUE
CALL KDRCPU(TK2)
IF(IMPX.GT.1) WRITE(6,'(/34H AUTFLU: CPU TIME SPENT TO COMPUTE,
1 31H SELF-SHIELDED REACTION RATES =,F8.1,19H SECOND, INCLUDING:
2 /9X,F8.1,36H SECOND TO SOLVE FOR AUTOSECOL FLUX;/9X,7HNUMBER ,
3 25HOF ASSEMBLY DOORS CALLS =,I5,1H.)') TK2-TK1,TK4,ICPIJ
*----
* SCRATCH STORAGE DEALLOCATION
*----
DEALLOCATE(III,NEXT,PRI)
DEALLOCATE(SCAT,SGAR,IJJ,NJJ)
DEALLOCATE(IPISO2,IPISO1)
DEALLOCATE(SIGINF,MASKH,DELTAU,CONC,GA2,GA1)
RETURN
540 FORMAT(/18H AUTFLU: ISOTOPE ',A12,1H'/12X,4HSIGS,16X,4HSIGT,
1 16X,4HSIGF,16X,5HSIGS1/(I5,1P,4E20.7))
550 FORMAT(/29H AUTFLU: SCATTERING ELEMENTS:,2I10/(1P,10E13.5))
END
|
