summaryrefslogtreecommitdiff
path: root/Dragon/src/SPHCMI.f
blob: d977cd3ce0f870097a8f42537defb9a1886332b4 (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
*DECK SPHCMI
      SUBROUTINE SPHCMI(IPMICR,IPRINT,IMC,NMERGE,NISOT,NGCOND,NL,NW,
     1 NED,NDEL,NALBP,SPH)
*
*-----------------------------------------------------------------------
*
*Purpose:
* SPH-correction of a Microlib.
*
*Copyright:
* Copyright (C) 2011 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
* IPMICR  pointer to the condensed microlib (L_LIBRARY signature).
* IPRINT  print flag (equal to 0 for no print).
* IMC     type of macro-calculation (=1 diffusion or SPN;
*         =2 other options;
*         =3 type PIJ with Bell acceleration).
* NMERGE  number of merged regions.
* NISOT   number of isotopes in microlib.
* NGCOND  number of condensed groups.
* NL      number of Legendre orders in scattering info.
* NW      type of weighting for PN cross section info (=0 P0; =1 P1).
* NED     number of additional phi-weighted edits in microlib.
* NDEL    number of delayed precursor groups.
* NALBP   number of physical albedos per condensed group.
* SPH     SPH homogenization factors.
*
*-----------------------------------------------------------------------
*
      USE GANLIB
*----
*  SUBROUTINE ARGUMENTS
*----
      TYPE(C_PTR) IPMICR
      INTEGER IPRINT,IMC,NMERGE,NISOT,NGCOND,NL,NW,NED,NDEL,NALBP
      REAL SPH(NMERGE+NALBP,NGCOND)
*----
*  LOCAL VARIABLES
*----
      PARAMETER(NSTATE=40)
      INTEGER ISTATE(NSTATE)
      TYPE(C_PTR) KPMICR
      CHARACTER HSIGN*12,TEXT12*12,CM*2,TEXT8*8,HSMG*131
      DOUBLE PRECISION DSUM
*----
*  ALLOCATABLE ARRAYS
*----
      INTEGER, ALLOCATABLE, DIMENSION(:) :: IMIX,ITYPR
      INTEGER, ALLOCATABLE, DIMENSION(:,:) :: IHUSED,IHEDIT
      REAL, ALLOCATABLE, DIMENSION(:,:) :: GAR
      REAL, ALLOCATABLE, DIMENSION(:,:,:) :: WSCAT
      CHARACTER(LEN=12), ALLOCATABLE, DIMENSION(:) :: HMAKE
      TYPE(C_PTR), ALLOCATABLE, DIMENSION(:) :: IPISO
*----
*  SCRATCH STORAGE ALLOCATION
*   IHEDIT  character*8 names of phi-weighted edits in microlib.
*----
      MAXH=7+2*NW+NL+NED+NDEL
      ALLOCATE(IHEDIT(2,NED+1),ITYPR(NL))
      ALLOCATE(GAR(NGCOND,MAXH),WSCAT(NGCOND,NGCOND,NL))
      ALLOCATE(IHUSED(3,NISOT),IMIX(NISOT),IPISO(NISOT))
      ALLOCATE(HMAKE(MAXH+NL))
*----
*  RECOVER MICROLIB INFORMATION
*----
      CALL LCMLEN(IPMICR,'SIGNATURE',LENGTH,ITYLCM)
      IF(LENGTH.EQ.0) GO TO 210
      CALL LCMGTC(IPMICR,'SIGNATURE',12,HSIGN)
      IF(HSIGN.NE.'L_LIBRARY') CALL XABORT('SPHCMI: MICROLIB EXPECTED')
      CALL LCMGET(IPMICR,'STATE-VECTOR',ISTATE)
      IF(ISTATE(1).NE.NMERGE) CALL XABORT('SPHCMI: INVALID NMERGE')
      IF(ISTATE(2).NE.NISOT) CALL XABORT('SPHCMI: INVALID NISOT')
      IF(ISTATE(3).NE.NGCOND) CALL XABORT('SPHCMI: INVALID NGCOND')
      IF(ISTATE(4).NE.NL) CALL XABORT('SPHCMI: INVALID NL')
      IF(ISTATE(13).NE.NED) CALL XABORT('SPHCMI: INVALID NED')
      IF(ISTATE(19).NE.NDEL) CALL XABORT('SPHCMI: INVALID NDEL')
      IF(NED.GT.0) CALL LCMGET(IPMICR,'ADDXSNAME-P0',IHEDIT)
*----
*  LOOP OVER ISOTOPES
*----
      CALL LCMGET(IPMICR,'ISOTOPESUSED',IHUSED)
      CALL LCMGET(IPMICR,'ISOTOPESMIX',IMIX)
      CALL LIBIPS(IPMICR,NISOT,IPISO)
      DO 200 ISOT=1,NISOT
      WRITE(TEXT12,'(3A4)') (IHUSED(I0,ISOT),I0=1,3)
      IF(IPRINT.GT.4) THEN
         WRITE(6,'(29H SPHCMI: PROCESSING ISOTOPE '',A12,2H''.)') TEXT12
      ENDIF
      IBM=IMIX(ISOT)
      KPMICR=IPISO(ISOT) ! set ISOT-th isotope
      IF(.NOT.C_ASSOCIATED(KPMICR)) THEN
         WRITE(HSMG,'(17HSPHCMI: ISOTOPE '',A12,16H'' IS NOT AVAILAB,
     >   19HLE IN THE MICROLIB.)') TEXT12
         CALL XABORT(HSMG)
      ENDIF
      MAXH=MAXH
      DO 10 J=1,MAXH+NL
      HMAKE(J)=' '
   10 CONTINUE
      GAR(:NGCOND,:MAXH)=0.0
      WSCAT(:NGCOND,:NGCOND,:NL)=0.0
*----
*  RECOVER CALCULATION-SPECIFIC ISOTOPIC DATA
*----
      DO 20 IW=1,MIN(NW+1,10)
      WRITE(TEXT12,'(3HNWT,I1)') IW-1
      CALL LCMLEN(KPMICR,TEXT12,LENGTH,ITYLCM)
      IF(LENGTH.EQ.NGCOND) THEN
         CALL LCMGET(KPMICR,TEXT12,GAR(1,IW))
      ELSE
         CALL LCMGET(KPMICR,'NWT0',GAR(1,IW))
      ENDIF
      HMAKE(IW)=TEXT12
      WRITE(TEXT12,'(4HNTOT,I1)') IW-1
      CALL LCMLEN(KPMICR,TEXT12,LENGTH,ITYLCM)
      IF(LENGTH.EQ.NGCOND) THEN
         CALL LCMGET(KPMICR,TEXT12,GAR(1,1+IW+NW))
      ELSE
         CALL LCMGET(KPMICR,'NTOT0',GAR(1,1+IW+NW))
      ENDIF
      HMAKE(1+IW+NW)=TEXT12
   20 CONTINUE
      CALL XDRLGS(KPMICR,-1,IPRINT,0,NL-1,1,NGCOND,GAR(1,3+2*NW),WSCAT,
     > ITYPR)
      DO 30 IL=0,NL-1
      IF(ITYPR(IL+1).NE.0) THEN
         WRITE (CM,'(I2.2)') IL
         HMAKE(3+2*NW+IL)='SIGS'//CM
      ENDIF
   30 CONTINUE
      CALL LCMLEN(KPMICR,'NUSIGF',LENGTH,ITYLCM)
      IF(LENGTH.EQ.NGCOND) THEN
         CALL LCMGET(KPMICR,'NUSIGF',GAR(1,3+2*NW+NL))
         HMAKE(3+2*NW+NL)='NUSIGF'
      ENDIF
      IF(NDEL.GT.0) THEN
         WRITE(TEXT12,'(6HNUSIGF,I2.2)') NDEL
         CALL LCMLEN(KPMICR,TEXT12,LENGTH,ITYLCM)
         IF(LENGTH.EQ.NGCOND) THEN
            DO 40 IDEL=1,NDEL
            WRITE(TEXT12,'(6HNUSIGF,I2.2)') IDEL
            CALL LCMGET(KPMICR,TEXT12,GAR(1,3+2*NW+NL+IDEL))
            HMAKE(3+2*NW+NL+IDEL)=TEXT12
   40       CONTINUE
         ENDIF
      ENDIF
      CALL LCMLEN(KPMICR,'H-FACTOR',LENGTH,ITYLCM)
      IF(LENGTH.EQ.NGCOND) THEN
         CALL LCMGET(KPMICR,'H-FACTOR',GAR(1,4+2*NW+NL+NDEL))
         HMAKE(4+2*NW+NL+NDEL)='H-FACTOR'
      ENDIF
      CALL LCMLEN(KPMICR,'OVERV',LENGTH,ITYLCM)
      IF(LENGTH.EQ.NGCOND) THEN
         CALL LCMGET(KPMICR,'OVERV',GAR(1,5+2*NW+NL+NDEL))
         HMAKE(5+2*NW+NL+NDEL)='OVERV'
      ENDIF
      CALL LCMLEN(KPMICR,'TRANC',LENGTH,ITYLCM)
      IF(LENGTH.EQ.NGCOND) THEN
         CALL LCMGET(KPMICR,'TRANC',GAR(1,6+2*NW+NL+NDEL))
         HMAKE(6+2*NW+NL+NDEL)='TRANC'
      ENDIF
      DO 60 IED=1,NED
      WRITE(TEXT8,'(2A4)') (IHEDIT(I0,IED),I0=1,2)
      CALL LCMLEN(KPMICR,TEXT8,LENGTH,ITYLCM)
      IF((LENGTH.GT.0).AND.(TEXT8.NE.'TRANC')) THEN
         CALL LCMGET(KPMICR,TEXT8,GAR(1,6+2*NW+NL+NDEL+IED))
         HMAKE(6+2*NW+NL+NDEL+IED)=TEXT8
      ENDIF
   60 CONTINUE
      CALL LCMLEN(KPMICR,'STRD',LENGTH,ITYLCM)
      IF(LENGTH.EQ.NGCOND) THEN
         CALL LCMGET(KPMICR,'STRD',GAR(1,7+2*NW+NL+NDEL+NED))
         HMAKE(7+2*NW+NL+NDEL+NED)='STRD'
      ENDIF
*----
*  APPLY SPH CORRECTION
*----
      DO 80 J=1,MAXH
      IF((HMAKE(J).NE.' ').AND.(HMAKE(J)(:4).NE.'SIGS')) THEN
         DO 70 IG=1,NGCOND
         IF((HMAKE(J)(:4).EQ.'NTOT').AND.(MOD(J-2-NW,2).EQ.1).AND.
     >      (IMC.EQ.1)) THEN
            GAR(IG,J)=GAR(IG,J)/SPH(IBM,IG)
         ELSE IF((HMAKE(J)(:3).EQ.'NWT').AND.(MOD(J-1,2).EQ.0)) THEN
            GAR(IG,J)=GAR(IG,J)/SPH(IBM,IG)
         ELSE IF((HMAKE(J)(:4).EQ.'STRD').OR.(HMAKE(J).EQ.'TRANC')) THEN
            GAR(IG,J)=GAR(IG,J)/SPH(IBM,IG)
         ELSE IF((HMAKE(J)(:3).EQ.'NWT').AND.(MOD(J-1,2).EQ.1)) THEN
            CONTINUE
         ELSE IF((HMAKE(J)(:4).EQ.'NTOT').AND.(MOD(J-2-NW,2).EQ.0).AND.
     >           (IMC.EQ.1)) THEN
            GAR(IG,J)=GAR(IG,J)*SPH(IBM,IG)
         ELSE IF((HMAKE(J)(:4).EQ.'NTOT').AND.(IMC.GT.1)) THEN
            CONTINUE
         ELSE
            GAR(IG,J)=GAR(IG,J)*SPH(IBM,IG)
         ENDIF
   70    CONTINUE
      ENDIF
   80 CONTINUE
      DO 105 IL=1,NL
      ITYPR(IL)=0
      IF(HMAKE(MAXH+IL+1).NE.' ') ITYPR(IL)=1
      DO 100 IG2=1,NGCOND
      DSUM=0.0
      DO 90 IG1=1,NGCOND
      IF(MOD(IL-1,2).EQ.0) THEN
         IF((IG1.EQ.IG2).AND.(IMC.GT.1).AND.(IL.LE.NW+1)) THEN
           WSCAT(IG1,IG1,IL)=WSCAT(IG1,IG1,IL)*SPH(IBM,IG1)
     >     +(GAR(IG1,2+NW)-GAR(IG1,1+IL+NW)*SPH(IBM,IG1))
         ELSE
           WSCAT(IG1,IG2,IL)=WSCAT(IG1,IG2,IL)*SPH(IBM,IG2) ! IG1 <- IG2
         ENDIF
      ELSE
         IF((IG1.EQ.IG2).AND.(IMC.GT.1).AND.(IL.LE.NW+1)) THEN
           WSCAT(IG1,IG1,IL)=WSCAT(IG1,IG1,IL)/SPH(IBM,IG1)
     >     +(GAR(IG1,2+NW)-GAR(IG1,1+IL+NW)/SPH(IBM,IG1))
         ELSE
           WSCAT(IG1,IG2,IL)=WSCAT(IG1,IG2,IL)/SPH(IBM,IG1)
         ENDIF
      ENDIF
      DSUM=DSUM+WSCAT(IG1,IG2,IL)
   90 CONTINUE
      IF((IL.EQ.1).AND.(IMC.GT.1)) THEN
         GAR(IG2,2+2*NW+IL)=GAR(IG2,2+2*NW+IL)*SPH(IBM,IG2)+
     >   GAR(IG2,2+NW)*(1.0-SPH(IBM,IG2))
      ELSE IF(IL.EQ.1) THEN
         GAR(IG2,2+2*NW+IL)=GAR(IG2,2+2*NW+IL)*SPH(IBM,IG2)
      ELSE
         GAR(IG2,2+2*NW+IL)=REAL(DSUM)
      ENDIF
  100 CONTINUE
  105 CONTINUE
*----
*  SAVE CORRECTED INFORMATION ON LCM
*----
      DO 110 J=1,MAXH
      IF((HMAKE(J).NE.' ').AND.(HMAKE(J)(:4).NE.'SIGS')) THEN
         CALL LCMPUT(KPMICR,HMAKE(J),NGCOND,2,GAR(1,J))
      ENDIF
  110 CONTINUE
      CALL XDRLGS(KPMICR,1,IPRINT,0,NL-1,1,NGCOND,GAR(1,3+2*NW),WSCAT,
     > ITYPR)
  200 CONTINUE
*----
*  CORRECT MACROLIB INFORMATION
*----
  210 CALL LCMLEN(IPMICR,'MACROLIB',LENGTH,ITYLCM)
      IF(LENGTH.NE.0) THEN
        CALL LCMSIX(IPMICR,'MACROLIB',1)
        CALL LCMLEN(IPMICR,'STATE-VECTOR',LENGTH,ITYLCM)
        IF(LENGTH.GT.0) THEN
          CALL LCMGET(IPMICR,'STATE-VECTOR',ISTATE)
          NIFISS=ISTATE(4)
          CALL SPHCMA(IPMICR,IPRINT,IMC,NMERGE,NGCOND,NIFISS,NED,NALBP,
     >    SPH)
        ENDIF
        CALL LCMSIX(IPMICR,' ',2)
      ENDIF
      IF(IPRINT.GT.5) WRITE(6,'(/28H SPHCMI: MICROLIB CORRECTED.)')
*----
*  SCRATCH STORAGE DEALLOCATION
*----
      DEALLOCATE(HMAKE)
      DEALLOCATE(IPISO,IMIX,IHUSED)
      DEALLOCATE(WSCAT,GAR)
      DEALLOCATE(ITYPR,IHEDIT)
      RETURN
      END