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
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
|
*DECK SHIDRV
SUBROUTINE SHIDRV (IPLIB,IPTRK,IFTRAK,LEVEL,NGRO,NBISO,NBMIX,
1 NREG,NUN,CDOOR,NRES,IMPX,ISONRF,ISONAM,MIX,DEN,SN,SB,LSHI,
2 IPHASE,IPROB,MAT,VOL,KEYFLX,LEAKSW,TITR,IGRMIN,IGRMAX,MAXX0,
3 IBIEFF,IGC,ITRANZ,EPS)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Perform a multidimensional self-shielding calculation in order to
* compute the dilution cross section of each resonant isotope present
* in the domain.
*
*Copyright:
* Copyright (C) 2002 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
* IPLIB pointer to the internal microscopic cross section library
* (L_LIBRARY signature).
* IPTRK pointer to the tracking. (L_TRACK signature).
* IFTRAK unit number of the sequential binary tracking file.
* LEVEL type of self-shielding model (=0 original Stamm'ler model;
* =1 original Stamm'ler model with Nordheim approximation;
* =2 Stamm'ler model with Nordheim approximation and Riemann
* integration method).
* NGRO number of energy groups.
* NBISO number of isotopes present in the calculation domain.
* NBMIX number of mixtures in the macrolib.
* NREG number of regions.
* NUN number of unknowns in the flux or source vector in one
* energy group.
* CDOOR name of the geometry/solution module.
* NRES number of resonant mixtures.
* IMPX print flag.
* ISONRF reference name of isotopes.
* ISONAM alias name of isotopes.
* MIX mix number of each isotope (can be zero).
* DEN density of each isotope.
* LSHI resonant region number associated with each isotope.
* Infinite dilution will be assumed if LSHI(i)=0.
* IPHASE type of flux solution (=1 use a native flux solution door;
* =2 use collision probabilities).
* IPROB adjoint macrolib flag (=0 direct; =1 adjoint).
* MAT index-number of the mixture type assigned to each volume.
* VOL volumes.
* KEYFLX pointers of fluxes in unknown vector.
* LEAKSW leakage flag (=.TRUE. only if leakage is present on the outer
* surface).
* TITR title.
* IGRMIN first group where the self-shielding is applied.
* IGRMAX most thermal group where the self-shielding is applied.
* MAXX0 maximum number of self-shielding iterations.
* IBIEFF Livolant-Jeanpierre normalization flag (=1 to activate).
* IGC Goldstein-Cohen approximation flag (=1 to activate).
* The Goldstein-Cohen approximation is activated only in cases
* where they are available in the internal library.
* ITRANZ type of transport correction used in the self-shielding
* calculations.
* EPS convergence criterion for the self-shielding iterations.
*
*Parameters: input/output
* SN on input, estimate of the dilution cross section in each
* energy group of each isotope. A value of 1.0e10 is used
* for infinite dilution and;
* on output, computed dilution cross section in each energy
* group of each isotope.
*
*Parameters: output
* SB dilution cross section as used in Livolant-Jeanpierre
* normalization.
*
*-----------------------------------------------------------------------
*
USE GANLIB
*----
* SUBROUTINE ARGUMENTS
*----
TYPE(C_PTR) IPLIB,IPTRK
INTEGER IFTRAK,LEVEL,NGRO,NBISO,NBMIX,NREG,NUN,NRES,IMPX,
1 ISONRF(3,NBISO),ISONAM(3,NBISO),MIX(NBISO),LSHI(NBISO),
2 IPHASE,IPROB,MAT(NREG),KEYFLX(NREG),IGRMIN,IGRMAX,MAXX0,IBIEFF,
3 IGC,ITRANZ
REAL DEN(NBISO),SN(NGRO,NBISO),SB(NGRO,NBISO),VOL(NREG),EPS
LOGICAL LEAKSW
CHARACTER CDOOR*12,TITR*72
*----
* LOCAL VARIABLES
*----
PARAMETER (NALPHA=9,NRAT=(NALPHA+1)/2,NSTATE=40)
TYPE(C_PTR) JPLIB,KPLIB
INTEGER IPAR(NSTATE)
REAL TMPDAY(3)
CHARACTER HSMG*130
LOGICAL START,BIEFF,LGC,LOGDO
*----
* ALLOCATABLE ARRAYS
*----
INTEGER, ALLOCATABLE, DIMENSION(:) :: LSHI2
REAL, ALLOCATABLE, DIMENSION(:) :: SIGE
REAL, ALLOCATABLE, DIMENSION(:,:) :: SIGT1,SIGT2,SIGT3
LOGICAL, ALLOCATABLE, DIMENSION(:) :: MASK,MASKL
LOGICAL, ALLOCATABLE, DIMENSION(:,:) :: NOCONV
*----
* SCRATCH STORAGE ALLOCATION
*----
ALLOCATE(SIGT1(NBMIX,NGRO),SIGT2(NBMIX,NGRO),SIGT3(NBMIX,NGRO))
ALLOCATE(MASK(NBMIX),MASKL(NGRO),NOCONV(NBMIX,NGRO))
*
IF(IMPX.GE.5) THEN
WRITE (6,'(//23H SHIDRV: VALUES OF MAT:)')
I1=1
KI=(NREG-1)/11+1
DO 10 I=1,KI
I2=I1+10
IF(I2.GT.NREG) I2=NREG
WRITE (6,340) (J,J=I1,I2)
WRITE (6,350) (MAT(J),J=I1,I2)
I1=I1+11
10 CONTINUE
WRITE (6,'(//)')
ENDIF
*----
* RECOVER SELF SHIELDING DATA
*----
IF(LEAKSW) CALL XABORT('SHIDRV: NEUTRON LEAKAGE IS FORBIDDEN.')
IF(CDOOR.EQ.' ') CALL XABORT('SHIDRV: THE GEOMETRY IS NOT YET '
1 //'DEFINED.')
BIEFF=(IBIEFF.EQ.1)
LGC=(IGC.EQ.1)
*
IF(IMPX.GT.0) THEN
WRITE (6,400) TITR,CDOOR
WRITE (6,'(25H STAMM''LER APPROXIMATION./)')
WRITE (6,405) IGRMIN,IGRMAX,MAXX0,IBIEFF,IGC,ITRANZ,LEVEL,
1 IPHASE,EPS
ENDIF
IF(NRES.EQ.0) THEN
WRITE (6,410)
RETURN
ENDIF
DO 30 I=1,NREG
IF(MAT(I).GT.NBMIX) THEN
WRITE (HSMG,380) NBMIX
CALL XABORT(HSMG)
ENDIF
30 CONTINUE
IGRMAX=MIN(IGRMAX,NGRO)
DO 55 LLL=1,NGRO
DO 40 I=1,NBISO
SB(LLL,I)=SN(LLL,I)
40 CONTINUE
DO 50 IBM=1,NBMIX
NOCONV(IBM,LLL)=.FALSE.
50 CONTINUE
55 CONTINUE
CALL LCMSIX(IPLIB,'MACROLIB',1)
JPLIB=LCMGID(IPLIB,'GROUP')
DO 70 LLL=IGRMIN,IGRMAX
IF(IPROB.EQ.0) LL=LLL
IF(IPROB.EQ.1) LL=NGRO-LLL+1
KPLIB=LCMGIL(JPLIB,LL)
CALL LCMGET(KPLIB,'NTOT0',SIGT2(1,LLL))
*----
* TRANSPORT CORRECTION
*----
IF(ITRANZ.NE.0) THEN
CALL LCMGET(KPLIB,'TRANC',SIGT3(1,LLL))
ELSE
SIGT3(:NBMIX,LLL)=0.0
ENDIF
*
NOCONV(:NBMIX,LLL)=.TRUE.
70 CONTINUE
CALL LCMSIX(IPLIB,' ',2)
IF(IMPX.GE.5) THEN
WRITE (6,'(/20H SHIBA INPUT VALUES:/)')
DO 80 L=IGRMIN,IGRMAX
WRITE(6,420) L
WRITE(6,460) (SN(L,J),J=1,NBISO)
WRITE(6,480) (SIGT2(IBM,L),IBM=1,NBMIX)
80 CONTINUE
WRITE(6,490)
ENDIF
*----
* ELIMINATE ISOTOPE ABSENT FROM GEOMETRY
*----
DO IBM=1,NBMIX
DO IREG=1,NREG
IF(MAT(IREG).EQ.IBM) GO TO 85
ENDDO
DO ISO=1,NBISO
IF(MIX(ISO).EQ.IBM) LSHI(ISO)=0
ENDDO
85 CONTINUE
ENDDO
*----
* RECOMPUTE THE VECTOR LSHI
*----
IF(LEVEL.EQ.0) THEN
NRES2=NRES
ELSE
ALLOCATE(LSHI2(NBISO))
NRES1=0
NRES2=0
DO 90 ISO=1,NBISO
LSHI2(ISO)=0
90 CONTINUE
DO 140 INRS=1,NRES
100 DENMAX=0.0
KSOT=0
DO 120 ISO=1,NBISO
IF(LSHI2(ISO).EQ.0) THEN
VOLISO=0.0
DO 110 I=1,NREG
IF(MAT(I).EQ.MIX(ISO)) VOLISO=VOLISO+VOL(I)
110 CONTINUE
IF((LSHI(ISO).EQ.INRS).AND.(DEN(ISO)*VOLISO.GT.DENMAX)) THEN
KSOT=ISO
DENMAX=DEN(ISO)*VOLISO
ENDIF
ENDIF
120 CONTINUE
IF(KSOT.GT.0) THEN
NRES2=NRES2+1
DO 130 ISO=1,NBISO
IF((ISONRF(1,ISO).EQ.ISONRF(1,KSOT)).AND.
1 (ISONRF(2,ISO).EQ.ISONRF(2,KSOT)).AND.
2 (ISONRF(3,ISO).EQ.ISONRF(3,KSOT)).AND.
3 (LSHI(ISO).EQ.INRS)) LSHI2(ISO)=NRES2
IF((ISONAM(1,ISO).EQ.ISONAM(1,KSOT)).AND.
1 (ISONAM(2,ISO).EQ.ISONAM(2,KSOT)).AND.
2 (LSHI(ISO).EQ.INRS)) LSHI2(ISO)=NRES2
130 CONTINUE
GO TO 100
ENDIF
IF(NRES2.EQ.NRES1) THEN
WRITE(HSMG,'(43HSHIDRV: NO RESONANT ISOTOPES IN RESONANT RE,
1 11HGION NUMBER,I4,5H (1).)') INRS
CALL XABORT(HSMG)
ENDIF
NRES1=NRES2
140 CONTINUE
ENDIF
*----
* DETERMINE THE AMOUNT OF SCRATCH STORAGE REQUIRED
*----
NBMIX2=0
DO 150 ISO=1,NBISO
IF(LSHI(ISO).GT.0) NBMIX2=NBMIX2+1
150 CONTINUE
*----
* ITERATION LOOP
*----
IF(LEVEL.EQ.0) ALLOCATE(SIGE(NRES2*NGRO))
NITER=0
160 NITER=NITER+1
START=(NITER.EQ.1)
IF(IMPX.GT.5) WRITE (6,430) NITER
DO 175 L=IGRMIN,IGRMAX
DO 170 IBM=1,NBMIX
SIGT1(IBM,L)=SIGT2(IBM,L)
170 CONTINUE
175 CONTINUE
IF(LEVEL.EQ.0) THEN
CALL SHISN2 (IPLIB,IPTRK,IFTRAK,NGRO,NBISO,NBMIX,NREG,NUN,
1 CDOOR,NRES,NBMIX2,IMPX,ISONAM,MIX,DEN,SN,SB,LSHI,IPHASE,
2 MAT,VOL,KEYFLX,LEAKSW,TITR,START,SIGT2,SIGT3,NOCONV,BIEFF,
3 LGC,SIGE)
ELSE
DO 210 INRS=1,NRES2
NBNRS=0
DO 200 IBM=1,NBMIX
LOGDO=.FALSE.
DO 180 I=1,NREG
LOGDO=LOGDO.OR.(MAT(I).EQ.IBM)
180 CONTINUE
IF(.NOT.LOGDO) GO TO 200
DO 190 ISO=1,NBISO
IF((MIX(ISO).EQ.IBM).AND.(LSHI2(ISO).EQ.INRS)) THEN
NBNRS=NBNRS+1
GO TO 200
ENDIF
190 CONTINUE
200 CONTINUE
IF(NBNRS.EQ.0) THEN
IF(START.AND.(IMPX.GE.1)) WRITE(6,385) 'SHIDRV',INRS
GO TO 210
ELSE IF(START.AND.(NBNRS.GT.1).AND.(IMPX.GE.1)) THEN
WRITE (6,370) NBNRS,INRS
ENDIF
CALL SHISN3 (IPLIB,IPTRK,IFTRAK,LEVEL,NGRO,NBISO,NBMIX,NREG,
1 NUN,CDOOR,INRS,NBNRS,IMPX,ISONAM,MIX,DEN,SN,SB,LSHI2,IPHASE,
2 MAT,VOL,KEYFLX,LEAKSW,TITR,START,SIGT2,SIGT3,NOCONV,BIEFF,LGC)
210 CONTINUE
ENDIF
ZZMAX=0.0
LNGRO=0
ICOUNT=0
DO 240 L=IGRMIN,IGRMAX
ZNORM=0.0
DO 220 IBM=1,NBMIX
ZNORM=MAX(ZNORM,ABS(SIGT2(IBM,L)))
220 CONTINUE
ZMAX=0.0
MASKL(L)=.FALSE.
DO 230 IBM=1,NBMIX
YMAX=ABS(SIGT1(IBM,L)-SIGT2(IBM,L))/ZNORM
ZMAX=MAX(ZMAX,YMAX)
NOCONV(IBM,L)=(NOCONV(IBM,L).AND.(YMAX.GT.EPS))
MASKL(L)=MASKL(L).OR.NOCONV(IBM,L)
230 CONTINUE
IF(MASKL(L)) ICOUNT=ICOUNT+1
IF(ZMAX.GT.ZZMAX) THEN
ZZMAX=ZMAX
LNGRO=L
ENDIF
240 CONTINUE
IF(IMPX.GE.2) WRITE (6,440) NITER,ICOUNT,ZZMAX,LNGRO
IF(IMPX.GE.10) THEN
WRITE (6,450) (L,MASKL(L),L=IGRMIN,IGRMAX)
WRITE (6,'(/31H INPUT MACROSCOPIC X-S IN GROUP,I4,1H:)') LNGRO
WRITE (6,'(1X,1P,10E12.4)') (SIGT1(IBM,LNGRO),IBM=1,NBMIX)
WRITE (6,'(/32H OUTPUT MACROSCOPIC X-S IN GROUP,I4,1H:)') LNGRO
WRITE (6,'(1X,1P,10E12.4)') (SIGT2(IBM,LNGRO),IBM=1,NBMIX)
ENDIF
IF(IMPX.GT.3) THEN
WRITE (6,'(/29H OUTPUT DILUTION X-S IN GROUP,I4,1H:)') LNGRO
WRITE (6,'(1X,1P,10E12.4)') (SN(LNGRO,J),J=1,NBISO)
ENDIF
IF((NITER.GE.MAXX0).AND.(ICOUNT.GT.0)) THEN
WRITE (6,390)
GO TO 250
ELSE IF(ICOUNT.GT.0) THEN
GO TO 160
ENDIF
*----
* CONVERGENCE IS OBTAINED
*----
250 IF(LEVEL.GT.0) DEALLOCATE(LSHI2)
IF(LEVEL.EQ.0) DEALLOCATE(SIGE)
IF(IMPX.GE.3) THEN
WRITE (6,'(/21H SHIBA OUTPUT VALUES:/)')
DO 260 L=IGRMIN,IGRMAX
WRITE(6,420) L
WRITE(6,460) (SN(L,J),J=1,NBISO)
IF(BIEFF) WRITE(6,470) (SB(L,J),J=1,NBISO)
IF(IMPX.GE.5) WRITE(6,480) (SIGT2(IBM,L),IBM=1,NBMIX)
260 CONTINUE
WRITE(6,490)
ENDIF
*----
* COMPUTE THE NEW SELF-SHIELDED MACROSCOPIC CROSS SECTIONS
*----
MASKL(:NGRO)=.FALSE.
DO 280 LLL=IGRMIN,IGRMAX
MASKL(LLL)=.TRUE.
280 CONTINUE
DO 300 IBM=1,NBMIX
DO 290 ISO=1,NBISO
MASK(IBM)=(MIX(ISO).EQ.IBM).AND.(LSHI(ISO).GT.0)
IF(MASK(IBM)) GO TO 300
290 CONTINUE
300 CONTINUE
ITSTMP=0
TMPDAY(1)=0.0
TMPDAY(2)=0.0
TMPDAY(3)=0.0
CALL LIBMIX(IPLIB,NBMIX,NGRO,NBISO,ISONAM,MIX,DEN,MASK,MASKL,
> ITSTMP,TMPDAY)
IF(IMPX.GT.0) WRITE (6,500) NITER,ZZMAX
*----
* STORE THE GENERAL SHIBA PARAMETERS
*----
IPAR(:NSTATE)=0
IPAR(1)=IGRMIN
IPAR(2)=IGRMAX
IPAR(3)=MAXX0
IPAR(4)=IBIEFF
IPAR(5)=IGC
IPAR(6)=ITRANZ
IPAR(7)=LEVEL
IPAR(8)=IPHASE
CALL LCMSIX(IPLIB,'SHIBA',1)
CALL LCMPUT(IPLIB,'STATE-VECTOR',NSTATE,1,IPAR)
CALL LCMPUT(IPLIB,'EPS-SHIBA',1,2,EPS)
CALL LCMSIX(IPLIB,' ',2)
*----
* SCRATCH STORAGE DEALLOCATION
*----
DEALLOCATE(NOCONV,MASKL,MASK)
DEALLOCATE(SIGT3,SIGT2,SIGT1)
RETURN
*
340 FORMAT(//26H VOLUME NB. :,11(I5,3X,1HI))
350 FORMAT( 26H MIXTURE (MAT) :,11(I5,3X,1HI))
370 FORMAT(/42H SHIDRV: USE THE NORDHEIM MODEL TO PROCESS,I3,5H RESO,
1 39HNANT MIXTURES IN RESONANT REGION NUMBER,I3,1H.)
380 FORMAT(32HSHIDRV: INVALID VALUE OF NBMIX (,I5,2H).)
385 FORMAT(A6,1X,': RESONANT REGION =',I10,1X,'NOT USED.')
390 FORMAT(/1X,61(1H*)/42H SHIDRV: MAXIMUM NUMBER OF SELF-SHIELDING ,
1 20HITERATIONS EXCEEDED./1X,61(1H*)/)
400 FORMAT(
> 1X,'SHIBA MULTIDIMENSIONAL SELF-SHIELDING CALCULATION',
> 1X,'-> A. HEBERT'/
> 1X,A72/
> 1X,'COLLISION PROBABILITY MODULE: ',A12/)
405 FORMAT(/8H OPTIONS/8H -------/
1 7H IGRMIN,I8,27H (FIRST GROUP TO PROCESS)/
2 7H IGRMAX,I8,34H (MOST THERMAL GROUP TO PROCESS)/
3 7H MAXX0 ,I8,33H (MAXIMUM NUMBER OF ITERATIONS)/
4 7H IBIEFF,I8,46H (=1: USE LIVOLANT-JEANPIERRE NORMALIZATION)/
5 7H IGC ,I8,42H (=1: USE GOLDSTEIN-COHEN APPROXIMATION)/
6 7H ITRANZ,I8,45H (0=NO TRANSPORT CORRECTION/1=APOLLO TYPE/2,
7 57H=RECOVER FROM LIBRARY/3=WIMS-D TYPE/4=LEAKAGE CORRECTION)/
8 7H LEVEL ,I8,46H (=0: STAMM'LER; =1: STAMM'LER/NORDHEIM; =2:,
9 18H RIEMANN/NORDHEIM)/
1 7H IPHASE,I8,37H (=1: NATIVE ASSEMBLY; =2: USE PIJ)/
2 7H EPS ,1P,E8.1,22H (STOPING CRITERION)/)
410 FORMAT(/52H SHIDRV: THERE IS NO REQUEST TO PROCESS ANY RESONANT,
1 9H ISOTOPE./)
420 FORMAT(/1X,131(1H-)//8H GROUP =,I4/)
430 FORMAT(/40H PERFORMING SELF-SHIELDING ITERATION NB.,I5)
440 FORMAT(/27H SELF-SHIELDING ITERATION =,I4,5X,14HNUMBER OF NON ,
1 18HCONVERGED GROUPS =,I4,5X,7HERROR =,1P,E13.4,0P,9H IN GROUP,
2 I4/)
450 FORMAT(7H MASKL(,I3,2H)=,L1,:,8H MASKL(,I3,2H)=,L1,:,8H MASKL(,
1 I3,2H)=,L1,:,8H MASKL(,I3,2H)=,L1,:,8H MASKL(,I3,2H)=,L1,:,
2 8H MASKL(,I3,2H)=,L1,:,8H MASKL(,I3,2H)=,L1,:,8H MASKL(,I3,
3 2H)=,L1,:,8H MASKL(,I3,2H)=,L1)
460 FORMAT(/37H MICROSCOPIC DILUTION CROSS SECTIONS:/(9X,1P,11E11.3))
470 FORMAT(/53H LIVOLANT AND JEANPIERRE MICROSCOPIC DILUTION CROSS S,
1 8HECTIONS:/(9X,1P,11E11.3))
480 FORMAT(/34H MACROSCOPIC TOTAL CROSS SECTIONS:/(9X,1P,11E11.3))
490 FORMAT(/1X,131(1H-)/)
500 FORMAT(/41H CONVERGENCE REACHED AT SHIBA ITERATION =,I4,7H ERROR,
1 2H =,1P,E11.3/)
END
|