summaryrefslogtreecommitdiff
path: root/Donjon/src/PCRMIC.f
blob: 9219e442fa05e93bf579d576b399f955432b2833 (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
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
*DECK PCRMIC
      SUBROUTINE PCRMIC(MAXNIS,MAXISO,IPLIB,IACCS,NMIX,NGRP,IMPX,
     1 NCAL,TERP,NISO,HISO,CONC,LMIXC,XS_CALC,B2)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Build the microlib by scanning the NCAL elementary calculations from
* PMAXS file and weighting them with TERP factors.
*
*Copyright:
* Copyright (C) 2019 Ecole Polytechnique de Montreal
*
*Author(s): 
* A. Hebert
*
*Parameters: input
* MAXNIS  maximum value of NISO(I) in user data.
* MAXISO  maximum allocated space for output microlib TOC information.
* IPLIB   address of the output microlib LCM object.
* IACCS   =0 microlib is created; =1 ... is updated.
* NMIX    maximum number of material mixtures in the microlib.
* NGRP    number of energy groups.
* IMPX    print parameter (equal to zero for no print).
* NCAL    number of elementary calculations in the PMAXS file.
* TERP    interpolation factors.
* NISO    number of user-selected isotopes.
* HISO    name of the user-selected isotopes.
* CONC    user-defined number density of the user-selected isotopes.
* LMIXC   flag set to .true. for fuel-map mixtures to process.
* XS_CALC pointers towards PMAXS elementary calculations.
* B2      buckling
*
*-----------------------------------------------------------------------
*
      USE GANLIB
      USE PCRDATA
      IMPLICIT NONE
*----
*  SUBROUTINE ARGUMENTS
*----
      TYPE(C_PTR) IPLIB
      INTEGER MAXNIS,MAXISO,IACCS,NMIX,NGRP,IMPX,NCAL,NISO(NMIX),
     1 HISO(2,NMIX,MAXNIS)
      REAL TERP(NCAL,NMIX),CONC(NMIX,MAXNIS),B2
      LOGICAL LMIXC(NMIX)
      TYPE(XSBLOCK_ITEM) XS_CALC(NCAL)
*----
*  LOCAL VARIABLES
*----
      INTEGER, PARAMETER::IOUT=6
      INTEGER, PARAMETER::MAXED=50
      INTEGER, PARAMETER::NSTATE=40
      INTEGER I0, IBM, ICAL, IED1, IED2, IGR, ISO, ITRANC, KSO1, I,
     & JSO, KSO, NBISO1, NBISO2, NCOMB2, NCOMB, NDEL, NDEPL, NED1,
     & NED2, NL, ITSTMP, MAXMIX, NBISO
      REAL WEIGHT,TMPDAY(3)
      CHARACTER TEXT12*12,HNAME*12,HVECT1(MAXED)*8,HHISO*8,TEXT8*8,
     & HVECT2(MAXED)*8
      INTEGER ISTATE(NSTATE)
      TYPE(C_PTR) IPTMP,JPTMP,KPTMP,JPLIB,KPLIB
*----
*  ALLOCATABLE ARRAYS
*----
      INTEGER, ALLOCATABLE, DIMENSION(:) :: IMIX2,ITOD2,ISTY2,MILVO,
     & IMICR
      INTEGER, ALLOCATABLE, DIMENSION(:,:) :: HUSE1,HUSE2,HNAM2
      REAL, ALLOCATABLE, DIMENSION(:) :: DENS2,ENER,GAR1,GAR2
      LOGICAL, ALLOCATABLE, DIMENSION(:) :: MASKL
      INTEGER, POINTER, DIMENSION(:) :: ISONA,ISOMI
      REAL, POINTER, DIMENSION(:) :: DENIS
      TYPE(C_PTR) ISONA_PTR,ISOMI_PTR,DENIS_PTR
      TYPE(C_PTR), ALLOCATABLE, DIMENSION(:) :: IPLIST
*----
*  SCRATCH STORAGE ALLOCATION
*----
      ALLOCATE(HUSE1(3,MAXISO),IMIX2(MAXISO),ITOD2(MAXISO),
     & ISTY2(MAXISO),HUSE2(3,MAXISO),HNAM2(3,MAXISO),MILVO(NMIX),
     & IPLIST(MAXISO))
      ALLOCATE(DENS2(MAXISO),ENER(NGRP+1))
*----
*  MICROLIB INITIALIZATION
*----
      ITRANC=0
      DENS2(:MAXISO)=0.0
      IMIX2(:MAXISO)=0
      ITOD2(:MAXISO)=0
      ISTY2(:MAXISO)=3
      IPLIST(:MAXISO)=C_NULL_PTR
      IF(IACCS.EQ.0) THEN
         NBISO2=0
         NCOMB2=0
         NED2=0
         TEXT12='L_LIBRARY'
         CALL LCMPTC(IPLIB,'SIGNATURE',12,TEXT12)
      ELSE
         CALL LCMGET(IPLIB,'STATE-VECTOR',ISTATE)
         IF(ISTATE(1).NE.NMIX) CALL XABORT('PCRMIC: INVALID NUMBER OF '
     1   //'MATERIAL MIXTURES IN THE MICROLIB.')
         IF(ISTATE(3).NE.NGRP) CALL XABORT('PCRMIC: INVALID NUMBER OF '
     1   //'ENERGY GROUPS IN THE MICROLIB.')
         NBISO2=ISTATE(2)
         NCOMB2=ISTATE(12)
         IF(NBISO2.GT.MAXISO) CALL XABORT('PCRMIC: MAXISO OVERFLOW(1).')
         NED2=ISTATE(13)
         IF(NED2.GT.MAXED) CALL XABORT('PCRMIC: MAXED OVERFLOW.')
         CALL LCMGET(IPLIB,'ISOTOPESUSED',HUSE2)
         CALL LCMGET(IPLIB,'ISOTOPERNAME',HNAM2)
         CALL LCMGET(IPLIB,'ISOTOPESDENS',DENS2)
         CALL LCMGET(IPLIB,'ISOTOPESMIX',IMIX2)
         CALL LCMGET(IPLIB,'ISOTOPESTODO',ITOD2)
         CALL LCMGET(IPLIB,'ISOTOPESTYPE',ISTY2)
         IF(NED2.GT.0) CALL LCMGTC(IPLIB,'ADDXSNAME-P0',8,NED2,HVECT2)
         CALL LCMGET(IPLIB,'ENERGY',ENER)
      ENDIF
*----
*  LOOP OVER MICROLIB MIXTURES
*----
      MILVO(:NMIX)=0
      NCOMB=0
      DO 190 IBM=1,NMIX
      IF(.NOT.LMIXC(IBM)) GO TO 190
      IF(NISO(IBM).GT.MAXNIS) CALL XABORT('PCRMIC: MAXNIS OVERFLOW.')
*----
*  FIND THE VALUE OF NBISO1 AND HUSE1 IN MIXTURE IBM
*----
      NBISO1=1
      TEXT12='*MAC*RES'
      READ(TEXT12,'(3A4)') (HUSE1(I,1),I=1,3)
      IF(NXST.GT.4) THEN
        NBISO1=3
        TEXT12='Xe135'
        READ(TEXT12,'(3A4)') (HUSE1(I,2),I=1,3)
        TEXT12='Sm149'
        READ(TEXT12,'(3A4)') (HUSE1(I,3),I=1,3)
      ENDIF
*----
*  LOOP OVER ELEMENTARY CALCULATIONS
*----
      CALL LCMOP(IPTMP,'*CALCULATIONS*',0,1,0)
      JPTMP=LCMLID(IPTMP,'CALCULATIONS',NCAL)
      DO 70 ICAL=1,NCAL
      WEIGHT=TERP(ICAL,IBM)
      IF(WEIGHT.EQ.0.0) GO TO 70
      KPTMP=LCMDIL(JPTMP,ICAL)
      CALL PCRONE(IMPX,ICAL,KPTMP,NCAL,NGRP,XS_CALC)
      IF(IMPX.GT.0) THEN
        WRITE(IOUT,'(33H PCRMIC: PMAXS ACCESS FOR MIXTURE,I8,6H AND C,
     1  10HALCULATION,I8,9H. WEIGHT=,1P,E12.4)') IBM,ICAL,WEIGHT
        IF(IMPX.GT.50) CALL LCMLIB(KPTMP)
      ENDIF
      CALL LCMGET(KPTMP,'STATE-VECTOR',ISTATE)
      IF(ISTATE(1).NE.1) CALL XABORT('PCRMIC: INVALID NUMBER OF MATERI'
     1 //'AL MIXTURES IN THE PMAXS FILE.')
      IF(ISTATE(2).NE.NBISO1) CALL XABORT('PCRMIC: INVALID NBISO1.')
      IF(ISTATE(3).NE.NGRP) CALL XABORT('PCRMIC: INVALID NUMBER OF ENE'
     1 //'RGY GROUPS IN THE COMPO.')
      NL=ISTATE(4)
      ITRANC=ISTATE(5)
      NDEPL=0
      NED1=ISTATE(13)
      NDEL=ISTATE(19)
      IF(NED1.GT.MAXED) CALL XABORT('PCRMIC: MAXED OVERFLOW.')
      IF(NED1.GT.0) CALL LCMGTC(KPTMP,'ADDXSNAME-P0',8,NED1,HVECT1)
      CALL LCMGET(KPTMP,'ENERGY',ENER)
      DO 30 IED1=1,NED1
      DO 20 IED2=1,NED2
      IF(HVECT1(IED1).EQ.HVECT2(IED2)) GO TO 30
   20 CONTINUE
      NED2=NED2+1
      HVECT2(NED2)=HVECT1(IED1)
   30 CONTINUE
      CALL LCMGPD(KPTMP,'ISOTOPESUSED',ISONA_PTR)
      CALL LCMGPD(KPTMP,'ISOTOPESDENS',DENIS_PTR)
      CALL C_F_POINTER(ISONA_PTR,ISONA,(/ NBISO1 /))
      CALL C_F_POINTER(DENIS_PTR,DENIS,(/ NBISO1 /))
      DO 60 ISO=1,NBISO1
      WRITE(TEXT8,'(2A4)') (ISONA(3*(ISO-1)+I0),I0=1,2)
      IF(TEXT8.EQ.'*MAC*RES') THEN
        DENIS(ISO)=1.0
      ELSE
        KSO1=0
        DO 40 KSO=1,NISO(IBM) ! user-selected isotope
        WRITE(HHISO,'(2A4)') (HISO(I0,IBM,KSO),I0=1,2)
        IF(TEXT8.EQ.HHISO) THEN
          KSO1=KSO
          GO TO 50
        ENDIF
   40   CONTINUE
   50   IF(KSO1.GT.0) DENIS(ISO)=CONC(IBM,KSO1)
      ENDIF
   60 CONTINUE
      CALL LCMPPD(KPTMP,'ISOTOPESDENS',NBISO1,2,DENIS_PTR)
   70 CONTINUE
*----
*  SELECT MICROLIB ISOTOPES CORRESPONDING TO PMAXS ISOTOPES
*----
      DO 90 ISO=1,NBISO1 ! PMAXS isotope
      WRITE(TEXT12,'(2A4)') (HUSE1(I,ISO),I=1,2)
      NBISO2=NBISO2+1
      IF(NBISO2.GT.MAXISO) THEN
        WRITE(IOUT,'(/16H PCRMIC: NBISO2=,I6,8H MAXISO=,I6)') NBISO2,
     1  MAXISO
        CALL XABORT('PCRMIC: MAXISO OVERFLOW(2).')
      ENDIF
      READ(TEXT12,'(3A4)') (HUSE2(I0,NBISO2),I0=1,3)
      DO 80 I0=1,3
      HNAM2(I0,NBISO2)=HUSE1(I0,ISO)
   80 CONTINUE
      IMIX2(NBISO2)=IBM
      DENS2(NBISO2)=DENIS(ISO)
   90 CONTINUE
      ALLOCATE(IMICR(NBISO1))
      IMICR(:NBISO1)=0
      DO 130 ISO=1,NBISO2 ! microlib isotope
      IF(IMIX2(ISO).NE.IBM) GO TO 130
      DO 120 JSO=1,NBISO1 ! PMAXS isotope
      IF((HUSE1(1,JSO).EQ.HUSE2(1,ISO)).AND.(HUSE1(2,JSO).EQ.
     1 HUSE2(2,ISO))) THEN
        IMICR(JSO)=ISO
        GO TO 130
      ENDIF
  120 CONTINUE
      WRITE(TEXT12,'(3A4)') (HUSE2(I0,ISO),I0=1,3)
      CALL XABORT('PCRMIC: UNABLE TO FIND '//TEXT12//'.')
  130 CONTINUE
*----
*  PROCESS ISOTOPE DIRECTORIES FOR MICROLIB MIXTURE IBM
*----
      DO 180 JSO=1,NBISO1 ! PMAXS isotope
      ISO=IMICR(JSO) ! microlib isotope
      IF(ISO.EQ.0) GO TO 180
      WRITE(HNAME,'(3A4)') (HUSE1(I0,JSO),I0=1,3)
      CALL LCMOP(KPLIB,'*ISOTOPE*',0,1,0)
      IPLIST(ISO)=KPLIB ! set isot ISO
      CALL PCRISO(KPLIB,JPTMP,HNAME,JSO,NCAL,NGRP,NL,NED2,HVECT2,NDEL,
     1 IMPX,TERP(1,IBM))
  180 CONTINUE
      DEALLOCATE(IMICR)
      CALL LCMCL(IPTMP,2)
  190 CONTINUE
*----
*  END OF LOOP OVER MICROLIB MIXTURES
*----
*----
*  CREATE ISOTOPE LIST DIRECTORY IN MICROLIB
*----
      JPLIB=LCMLID(IPLIB,'ISOTOPESLIST',NBISO2)
      DO 195 ISO=1,NBISO2 ! microlib isotope
      IF(C_ASSOCIATED(IPLIST(ISO))) THEN
        KPLIB=LCMDIL(JPLIB,ISO) ! step up isot ISO
        CALL LCMEQU(IPLIST(ISO),KPLIB)
        CALL LCMCL(IPLIST(ISO),2)
        IPLIST(ISO)=C_NULL_PTR
      ENDIF
  195 CONTINUE
*----
*  MICROLIB FINALIZATION
*----
      ISTATE(:NSTATE)=0
      ISTATE(1)=NMIX
      ISTATE(2)=NBISO2
      ISTATE(3)=NGRP
      ISTATE(4)=NL
      ISTATE(5)=ITRANC
      ISTATE(7)=1
      ISTATE(11)=NDEPL
      ISTATE(12)=NCOMB+NCOMB2
      ISTATE(13)=NED2
      ISTATE(14)=NMIX
      ISTATE(18)=1
      ISTATE(19)=NDEL
      ISTATE(22)=MAXISO/NMIX
      IF(NBISO2.EQ.0) CALL XABORT('PCRMIC: NBISO2=0.')
      CALL LCMPUT(IPLIB,'STATE-VECTOR',NSTATE,1,ISTATE)
      CALL LCMPUT(IPLIB,'ISOTOPESUSED',3*NBISO2,3,HUSE2)
      CALL LCMPUT(IPLIB,'ISOTOPERNAME',3*NBISO2,3,HNAM2)
      CALL LCMPUT(IPLIB,'ISOTOPESDENS',NBISO2,2,DENS2)
      CALL LCMPUT(IPLIB,'ISOTOPESMIX',NBISO2,1,IMIX2)
      IF(NED2.GT.0) CALL LCMPTC(IPLIB,'ADDXSNAME-P0',8,NED2,HVECT2)
      CALL LCMPUT(IPLIB,'ISOTOPESTODO',NBISO2,1,ITOD2)
      CALL LCMPUT(IPLIB,'ISOTOPESTYPE',NBISO2,1,ISTY2)
      CALL LCMPUT(IPLIB,'ENERGY',NGRP+1,2,ENER)
      IF(IMPX.GT.5) CALL LCMLIB(IPLIB)
      IACCS=1
*----
*  BUILD EMBEDDED MACROLIB
*----
      CALL LCMGET(IPLIB,'STATE-VECTOR',ISTATE)
      MAXMIX=ISTATE(1)
      IF(MAXMIX.NE.NMIX) CALL XABORT('PCRMIC: INVALID NMIX.')
      NBISO=ISTATE(2)
      ALLOCATE(MASKL(NGRP))
      CALL LCMGPD(IPLIB,'ISOTOPESUSED',ISONA_PTR)
      CALL LCMGPD(IPLIB,'ISOTOPESMIX',ISOMI_PTR)
      CALL LCMGPD(IPLIB,'ISOTOPESDENS',DENIS_PTR)
      CALL C_F_POINTER(ISONA_PTR,ISONA,(/ NBISO /))
      CALL C_F_POINTER(ISOMI_PTR,ISOMI,(/ NBISO /))
      CALL C_F_POINTER(DENIS_PTR,DENIS,(/ NBISO /))
      MASKL(:NGRP)=.TRUE.
      ITSTMP=0
      TMPDAY(1)=0.0
      TMPDAY(2)=0.0
      TMPDAY(3)=0.0
      CALL LIBMIX(IPLIB,MAXMIX,NGRP,NBISO,ISONA,ISOMI,DENIS,LMIXC,MASKL,
     1 ITSTMP,TMPDAY)
      DEALLOCATE(MASKL)
*----
*  INCLUDE LEAKAGE IN THE MACROLIB (USED ONLY FOR NON-REGRESSION TESTS)
*----
      IF(B2.NE.0.0) THEN
        IF(IMPX.GT.0) WRITE(6,'(/34H PCRMIC: INCLUDE LEAKAGE IN THE MA,
     1  11HCROLIB (B2=,1P,E12.5,2H).)') B2
        CALL LCMSIX(IPLIB,'MACROLIB',1)
        JPLIB=LCMGID(IPLIB,'GROUP')
        ALLOCATE(GAR1(NMIX),GAR2(NMIX))
        DO 210 IGR=1,NGRP
          KPLIB=LCMGIL(JPLIB,IGR)
          CALL LCMGET(KPLIB,'NTOT0',GAR1)
          CALL LCMGET(KPLIB,'DIFF',GAR2)
          DO 200 IBM=1,NMIX
            IF(LMIXC(IBM)) GAR1(IBM)=GAR1(IBM)+B2*GAR2(IBM)
  200     CONTINUE
          CALL LCMPUT(KPLIB,'NTOT0',NMIX,2,GAR1)
  210   CONTINUE
        DEALLOCATE(GAR2,GAR1)
        CALL LCMSIX(IPLIB,' ',2)
      ENDIF
*----
*  SCRATCH STORAGE DEALLOCATION
*----
      DEALLOCATE(ENER,DENS2)
      DEALLOCATE(IPLIST,MILVO,HNAM2,HUSE2,ISTY2,ITOD2,IMIX2,HUSE1)
      RETURN
      END