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
|
*DECK SPHDRV
SUBROUTINE SPHDRV(IPTRK,IFTRK,IPMACR,IPFLX,IPRINT,IMC,NGCOND,
1 NMERGE,NALBP,IGRMIN,IGRMAX,SPH)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Calculation of the SPH factors. These factors are used to multiply the
* cross sections and to divide the fluxes. The SPH factors calculation
* is generally application dependent. New SPH algorithms should be
* implemented in this driver.
*
*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
* IPTRK pointer to the macro-tracking LCM object.
* IFTRK unit of the macro-tracking binary sequential file.
* IPMACR pointer to the Macrolib (L_MACROLIB signature).
* IPFLX pointer towards an initialization flux (L_FLUX 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).
* NGCOND number of condensed groups.
* NMERGE number of merged regions.
* NALBP number of physical albedos.
* IGRMIN first group to process.
* IGRMAX last group to process.
*
*Parameters: output
* SPH SPH homogenization factors.
*
*-----------------------------------------------------------------------
*
USE GANLIB
*----
* SUBROUTINE ARGUMENTS
*----
TYPE(C_PTR) IPTRK,IPMACR,IPFLX
INTEGER IFTRK,IPRINT,IMC,NGCOND,NMERGE,NALBP,IGRMIN,IGRMAX
REAL SPH(NMERGE+NALBP,NGCOND)
*----
* LOCAL VARIABLES
*----
PARAMETER (NSTATE=40)
CHARACTER TEXT12*12,CNDOOR*12,CTITRE*72,SUFF(2)*2,HSMG*131
INTEGER ISTATE(NSTATE)
LOGICAL ILK
TYPE(C_PTR) IPSPH
*----
* ALLOCATABLE ARRAYS
*----
INTEGER, ALLOCATABLE, DIMENSION(:) :: MAT2,KEY2,MERG2
REAL, ALLOCATABLE, DIMENSION(:) :: VOL2,VREF,VMAC
DATA SUFF/'00','01'/
*----
* RECOVER SPH-RELATED INFORMATION
*----
CALL LCMLEN(IPMACR,'SPH',ILONG,ITYLCM)
IF(ILONG.EQ.0) CALL XABORT('SPHDRV: MISSING SPH DIRECTORY.')
IPSPH=LCMDID(IPMACR,'SPH')
CALL LCMGET(IPSPH,'STATE-VECTOR',ISTATE)
NSPH=ISTATE(1)
KSPH=ISTATE(2)
MAXIT=ISTATE(3)
MAXNBI=ISTATE(4)
IF((NSPH.EQ.0).OR.(NSPH.EQ.1)) CALL XABORT('SPHDRV: INVALID VALU'
> //'E OF NSPH.')
*----
* RECOVER AND USE AN EXISTING MACRO-TRACKING.
*----
IF(C_ASSOCIATED(IPTRK)) THEN
IF(NSPH.GE.2) THEN
CALL LCMGTC(IPSPH,'SPH$TRK',12,CNDOOR)
CALL LCMGET(IPSPH,'SPH-EPSILON',EPSPH)
ENDIF
CALL LCMGTC(IPTRK,'SIGNATURE',12,TEXT12)
IF(TEXT12.NE.'L_TRACK') THEN
CALL XABORT('SPHDRV: TRACKING DATA STRUCTURE EXPECTED.')
ENDIF
CALL LCMGET(IPTRK,'STATE-VECTOR',ISTATE)
NREG2=ISTATE(1)
NUN2=ISTATE(2)
ILK=ISTATE(3).EQ.0
CALL LCMGTC(IPTRK,'TRACK-TYPE',12,CNDOOR)
IF(CNDOOR.EQ.'MCCG') THEN
CALL LCMLEN(IPTRK,'KEYFLX',LKFL,ITYLCM)
NFUNL=LKFL/NREG2
ELSE
NFUNL=1
ENDIF
*----
* THE NUMBER OF UNKNOWNS OF A CURRENT-BASED MULTICELL ITERATION IS
* INCREASED TO HOLD INTERFACE CURRENT COMPONENTS.
*----
IF(NSPH.EQ.4) THEN
IF(CNDOOR.EQ.'SYBIL') NUN2=NUN2+ISTATE(9)
IF((CNDOOR.EQ.'EXCELL').AND.(ISTATE(7).EQ.5))
> NUN2=NUN2+ISTATE(28)
ENDIF
*
IF((CNDOOR.EQ.'EXCELL').OR.(CNDOOR.EQ.'MCCG')) THEN
ISCAT=ISTATE(6)
ELSE IF(CNDOOR.EQ.'SN') THEN
ISCAT=ISTATE(16)
ELSE IF(CNDOOR.EQ.'BIVAC') THEN
IF(ISTATE(14).EQ.0) THEN
ISCAT=1
ELSE
ISCAT=ISTATE(16)
ENDIF
ELSE IF(CNDOOR.EQ.'TRIVAC') THEN
IF(ISTATE(30).EQ.0) THEN
ISCAT=1
ELSE
ISCAT=ISTATE(32)
ENDIF
ELSE
ISCAT=1
ENDIF
ISCAT=ABS(ISCAT)
ALLOCATE(VOL2(NREG2),MAT2(NREG2),KEY2(NREG2*NFUNL))
CALL LCMGET(IPTRK,'VOLUME',VOL2)
CALL LCMGET(IPTRK,'MATCOD',MAT2)
CALL LCMGET(IPTRK,'KEYFLX',KEY2)
CALL LCMLEN(IPTRK,'TITLE',LENGT,ITYLCM)
IF(LENGT.GT.0) THEN
CALL LCMGTC(IPTRK,'TITLE',72,CTITRE)
ELSE
CTITRE='*** NO TITLE PROVIDED ***'
ENDIF
NBMIX2=0
IF(KSPH.EQ.5) THEN
* HEBERT-BENOIST ALBS TECHNIQUE.
DO 20 IREG=1,NREG2
NBMIX2=MAX(NBMIX2,MAT2(IREG))
20 CONTINUE
ALLOCATE(MERG2(NBMIX2))
DO 25 IBM=1,NBMIX2
MERG2(IBM)=1
25 CONTINUE
ILK=.FALSE.
ELSE
DO 30 IREG=1,NREG2
NBMIX2=MAX(NBMIX2,MAT2(IREG))
30 CONTINUE
IF(NBMIX2.NE.NMERGE) THEN
WRITE(HSMG,'(41HSPHDRV: INVALID NUMBER OF MACRO-REGIONS (,
> 2I6,2H).)') NBMIX2,NMERGE
CALL XABORT(HSMG)
ENDIF
ALLOCATE(MERG2(NBMIX2))
DO 35 IBM=1,NBMIX2
MERG2(IBM)=IBM
35 CONTINUE
ENDIF
*
* RECOVER TABULATED FUNCTIONS.
CALL XDRTA2
ELSE
ISCAT=1
ILK=.FALSE.
NBMIX2=0
NREG2=0
NUN2=0
ENDIF
*----
* CHECK VOLUME CONSISTENCY
*----
ALLOCATE(VREF(NMERGE),VMAC(NMERGE))
VMAC(:)=0.0
CALL LCMGET(IPMACR,'VOLUME',VREF)
DO IREG=1,NREG2
IBM=MAT2(IREG)
IF(IBM.GT.0) VMAC(IBM)=VMAC(IBM)+VOL2(IREG)
ENDDO
VREFT=SUM(VREF(:NMERGE))
VMACT=SUM(VMAC(:NMERGE))
DO IBM=1,NMERGE
ERR=ABS(VREF(IBM)/VREFT-VMAC(IBM)/VMACT)
IF(ERR.GT.1.0E-4*ABS(VREF(IBM)/VREFT)) THEN
WRITE(HSMG,'(38HSPHDRV: INCONSISTENT VOLUME IN MIXTURE,I5,
> 3H BY,F7.2,2H %)') IBM,ERR*100.0
CALL XABORT(HSMG)
ENDIF
ENDDO
DEALLOCATE(VMAC,VREF)
*----
* GENERAL PROCEDURE FOR COMPUTING THE SPH FACTORS
*----
CALL SPHEQU(NBMIX2,IPTRK,IFTRK,IPMACR,IPFLX,CNDOOR,NSPH,KSPH,
1 MAXIT,MAXNBI,EPSPH,IPRINT,IMC,NGCOND,NMERGE,NALBP,ISCAT,NREG2,
2 NUN2,MAT2,VOL2,KEY2,MERG2,ILK,CTITRE,IGRMIN,IGRMAX,SPH)
IF(C_ASSOCIATED(IPTRK)) DEALLOCATE(MERG2,KEY2,MAT2,VOL2)
*----
* PRINT SPH FACTORS
*----
IF(IPRINT.GT.1) THEN
WRITE(6,'(/21H SPHDRV: SPH FACTORS:)')
WRITE(6,200) ((IKK,IGR,SPH(IKK,IGR),IKK=1,NMERGE+NALBP),IGR=1,
> NGCOND)
ENDIF
RETURN
*
200 FORMAT(4X,4HSPH(,I5,1H,,I3,2H)=,F9.5,:,4X,4HSPH(,I5,1H,,I3,2H)=,
> F9.5,:,4X,4HSPH(,I5,1H,,I3,2H)=,F9.5,:,4X,4HSPH(,I5,1H,,I3,2H)=,
> F9.5,:,4X,4HSPH(,I5,1H,,I3,2H)=,F9.5)
END
|