summaryrefslogtreecommitdiff
path: root/Dragon/src/TONSN3.f
blob: e5ced06e0a4c19cc8f14d98e6e8a90818682ee14 (plain)
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
*DECK TONSN3
      SUBROUTINE TONSN3 (IPLIB,IPTRK,IFTRAK,NGRO,NBISO,NBM,NREG,NUN,
     1 CDOOR,INRS,NBNRS,IMPX,ISONAM,MIX,DEN,SN,LSHI,IPHASE,MAT,VOL,
     2 KEYFLX,LEAKSW,TITR,START,SIGT2,SIGT3,NOCONV,ICPIJ,TK3,TK4)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Perform one multidimensional self-shielding iteration using the
* Tone's method with Nordheim (PIC) approximation.
*
*Copyright:
* Copyright (C) 2017 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.
* NGRO    number of energy groups.
* NBISO   number of isotopes present in the calculation domain.
* NBM     number of mixtures in the macrolib.
* NREG    number of volumes.
* NUN     number of unknowns in the flux or source vector in one
*         energy group.
* CDOOR   name of the geometry/solution module.
* INRS    index of the resonant isotope under consideration.
* NBNRS   number of totaly correlated resonant regions.
* IMPX    print flag.
* 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).
* MAT     index-number of the mixture type assigned to each volume.
* VOL     volumes.
* KEYFLX  pointers of fluxes in unknown vector.
* LEAKSW  leakage flag (=.TRUE. if leakage is present on the outer
*         surface).
* TITR    title.
* START   beginning-of-iteration flag (=.TRUE. if TONSN3 is called
*         for the first time).
* SIGT3   transport correction.
* NOCONV  mixture convergence flag. (NOCONV(IBM,L)=.TRUE. if mixture IBM
*         is not converged in group L).
*
*Parameters: input/output
* SN      estimate of the dilution cross section in each energy group
*         of each isotope on input and computed dilution cross section
*         in each energy group of each isotope at output.
* SIGT2   total macroscopic cross sections on ipput and total 
*         macroscopic cross sections as modified by Tone's method
*         at output.
*
*Parameters: output
* ICPIJ   number of flux solution door calls.
*
*Reference:
* A. Hebert, 'Revisiting the Stamm'ler Self-Shielding Method', Presented
* at the 25th CNS annnual conference, June 6-9, Toronto, 2004.
*
*-----------------------------------------------------------------------
*
      USE GANLIB
*----
*  SUBROUTINE ARGUMENTS
*----
      TYPE(C_PTR) IPLIB,IPTRK
      INTEGER IFTRAK,NGRO,NBISO,NBM,NREG,NUN,INRS,NBNRS,IMPX,
     1 ISONAM(3,NBISO),MIX(NBISO),LSHI(NBISO),IPHASE,MAT(NREG),
     2 KEYFLX(NREG),ICPIJ
      REAL DEN(NBISO),SN(NGRO,NBISO),VOL(NREG),SIGT2(0:NBM,NGRO),
     1 SIGT3(0:NBM,NGRO),TK3,TK4
      CHARACTER CDOOR*12,TITR*72
      LOGICAL LEAKSW,START,NOCONV(NBM,NGRO)
*----
*  LOCAL VARIABLES
*----
      TYPE(C_PTR) KPLIB
      CHARACTER TEXT12*12,HNAMIS*12
      LOGICAL LOGDO
*----
*  ALLOCATABLE ARRAYS
*----
      INTEGER, ALLOCATABLE, DIMENSION(:) :: IRES,ISONR,NPSYS
      REAL, ALLOCATABLE, DIMENSION(:) :: GAR,GAS,VST,DENM
      REAL, ALLOCATABLE, DIMENSION(:,:) :: SIGT0,TOTAL,SIGE
      LOGICAL, ALLOCATABLE, DIMENSION(:) :: MASKI
      TYPE(C_PTR), ALLOCATABLE, DIMENSION(:) :: IPISO
*----
*  SCRATCH STORAGE ALLOCATION
*   SIGT0   macroscopic xs of the resonant isotopes as interpolated.
*----
      ALLOCATE(IRES(NBM),ISONR(NBISO),NPSYS(NGRO))
      ALLOCATE(SIGT0(0:NBM,NGRO),TOTAL(NGRO,NBNRS),DENM(0:NBM),
     1 GAR(NGRO),GAS(NGRO),SIGE(NBNRS,NGRO),VST(NBNRS))
      ALLOCATE(MASKI(NBISO))
      ALLOCATE(IPISO(NBISO))
*----
*  FIND THE RESONANT MIXTURE NUMBERS AND THE CORRELATED ISOTOPES
*  ASSOCIATED WITH REGION INRS
*----
      IRES(:NBM)=0
      ISONR(:NBISO)=0
      MASKI(:NBISO)=.FALSE.
      IRS=0
      TEXT12=' '
      DO 30 IBM=1,NBM
      LOGDO=.FALSE.
      DO 10 I=1,NREG
      LOGDO=LOGDO.OR.(MAT(I).EQ.IBM)
   10 CONTINUE
      IF(.NOT.LOGDO) GO TO 30
      DO 20 ISO=1,NBISO
      LOGDO=START.OR.(DEN(ISO).NE.0.)
      IF((MIX(ISO).EQ.IBM).AND.(LSHI(ISO).EQ.INRS)) THEN
         WRITE(HNAMIS,'(3A4)') (ISONAM(I0,ISO),I0=1,3)
         IF(HNAMIS.NE.TEXT12) THEN
           IRS=IRS+1
           TEXT12=HNAMIS
           IF(LOGDO) MASKI(ISO)=.TRUE.
         ENDIF
         ISONR(ISO)=IRS
         IRES(IBM)=IRS
      ENDIF
   20 CONTINUE
   30 CONTINUE
      IF(IRS.NE.NBNRS) CALL XABORT('TONSN3: INVALID VALUE OF NBNRS.')
*----
*  SET THE LCM MICROLIB ISOTOPEWISE DIRECTORIES.
*----
      CALL LIBIPS(IPLIB,NBISO,IPISO)
*----
*  UNLOAD MICROSCOPIC X-S FROM LCM TO SCRATCH STORAGE.
*----
      DO 40 ISO=1,NBISO
      IRS=ISONR(ISO)
      IF(IRS.GT.0) THEN
        KPLIB=IPISO(ISO) ! set ISO-th isotope
        CALL LCMGET(KPLIB,'NTOT0',TOTAL(1,IRS))
      ENDIF
   40 CONTINUE
*
      VST(:NBNRS)=0.0
      DO 60 I=1,NREG
      IF(MAT(I).EQ.0) GO TO 60
      IND=IRES(MAT(I))
      IF(IND.GT.0) VST(IND)=VST(IND)+VOL(I)
   60 CONTINUE
*
      NPSYS(:NGRO)=0
      DO 110 LLL=1,NGRO
      LOGDO=.FALSE.
      DO 70 IBM=1,NBM
      IRS=IRES(IBM)
      IF(IRS.GT.0) LOGDO=LOGDO.OR.NOCONV(IBM,LLL)
   70 CONTINUE
      IF(LOGDO) THEN
         NPSYS(LLL)=LLL
*
*        COMPUTE THE LIGHT AND RESONANT COMPONENTS OF THE MACROSCOPIC
*        CROSS SECTIONS IN EACH RESONANT MIXTURE.
         DENM(0:NBM)=0.0
         SIGT0(0:NBM,LLL)=0.0
         DO 90 ISO=1,NBISO
         IRS=ISONR(ISO)
         IF(IRS.GT.0) THEN
           IBM=MIX(ISO)
           DENM(IBM)=DEN(ISO)
           SIGT2(IBM,LLL)=SIGT2(IBM,LLL)-TOTAL(LLL,IRS)*DEN(ISO)
           SIGT0(IBM,LLL)=TOTAL(LLL,IRS)*DEN(ISO)
         ENDIF
   90    CONTINUE
         IF(IMPX.GE.10) THEN
            WRITE (6,400) LLL,(SIGT0(I,LLL),I=1,NBM)
            WRITE (6,410) LLL,(SIGT2(I,LLL),I=1,NBM)
         ENDIF
      ENDIF
  110 CONTINUE
*----
*  SET UP VECTOR SIGE.
*----
      CALL LCMSIX(IPLIB,'SHIBA',1)
*
      SIGE(:NBNRS,:NGRO)=0.0
      CALL LCMSIX(IPLIB,'--AVERAGE--',1)
      CALL TONDST(IPLIB,NPSYS,IPTRK,IFTRAK,CDOOR,IMPX,NBM,NBNRS,NREG,
     1 NUN,NGRO,IPHASE,MAT,VOL,KEYFLX,LEAKSW,IRES,DENM,SIGT0,SIGT2,
     2 SIGT3,TITR,SIGE,TK3,TK4)
      CALL LCMSIX(IPLIB,' ',2)
      DO 130 LLL=1,NGRO
      IF(NPSYS(LLL).NE.0) THEN
         ICPIJ=ICPIJ+2
         DO 120 ISO=1,NBISO
         IRS=ISONR(ISO)
         IF((LSHI(ISO).EQ.INRS).AND.(IRS.GT.0).AND.
     1   (DEN(ISO).NE.0.0)) THEN
            SN(LLL,ISO)=MAX(1.0,SIGE(IRS,LLL))
         ELSE IF((LSHI(ISO).EQ.INRS).AND.(IRS.GT.0)) THEN
            SN(LLL,ISO)=1.0E10
         ENDIF
  120    CONTINUE
      ENDIF
  130 CONTINUE
      CALL LCMSIX(IPLIB,' ',2)
*
      DO 320 LLL=1,NGRO
      LOGDO=.FALSE.
      DO 300 IBM=1,NBM
      IRS=IRES(IBM)
      IF(IRS.GT.0) LOGDO=LOGDO.OR.NOCONV(IBM,LLL)
  300 CONTINUE
      IF(LOGDO) THEN
         DO 310 ISO=1,NBISO
         IRS=ISONR(ISO)
         IF(IRS.GT.0) THEN
           IBM=MIX(ISO)
           SIGT2(IBM,LLL)=SIGT2(IBM,LLL)+TOTAL(LLL,IRS)*DEN(ISO)
         ENDIF
  310    CONTINUE
      ENDIF
  320 CONTINUE
*----
*  SAVE THE GROUP- AND ISOTOPE-DEPENDENT DILUTIONS
*----
      CALL LCMPUT(IPLIB,'ISOTOPESDSB',NBISO*NGRO,2,SN)
      CALL LCMPUT(IPLIB,'ISOTOPESDSN',NBISO*NGRO,2,SN)
*----
*  COMPUTE THE SELF-SHIELDED MICROSCOPIC CROSS SECTIONS AND UPDATE
*  VECTOR SIGT2
*----
      IMPX2=MAX(0,IMPX-1)
      CALL LIBLIB (IPLIB,NBISO,MASKI,IMPX2)
      DO 350 ISO=1,NBISO
      IRS=ISONR(ISO)
      IF(IRS.GT.0) THEN
        IBM=MIX(ISO)
        KPLIB=IPISO(ISO) ! set ISO-th isotope
        CALL LCMGET(KPLIB,'NTOT0',GAR)
        DO 340 LLL=1,NGRO
        TOTAL(LLL,IRS)=TOTAL(LLL,IRS)-GAR(LLL)
        SIGT2(IBM,LLL)=SIGT2(IBM,LLL)-DEN(ISO)*TOTAL(LLL,IRS)
  340   CONTINUE
      ENDIF
  350 CONTINUE
*----
*  SCRATCH STORAGE DEALLOCATION
*----
      DEALLOCATE(IPISO)
      DEALLOCATE(MASKI)
      DEALLOCATE(VST,SIGE,GAS,GAR,DENM,TOTAL,SIGT0)
      DEALLOCATE(NPSYS,ISONR,IRES)
      RETURN
*
  400 FORMAT(//51H TOTAL MACROSCOPIC CROSS SECTIONS OF THE RESONANT M,
     1 31HATERIALS IN EACH MIXTURE (GROUP,I5,2H):/(1X,1P,11E11.3))
  410 FORMAT(//51H TOTAL MACROSCOPIC CROSS SECTIONS OF THE OTHER MATE,
     1 28HRIALS IN EACH MIXTURE (GROUP,I5,2H):/(1X,1P,11E11.3))
      END