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
|
*DECK LIBTE2
SUBROUTINE LIBTE2 (NGRO,NSUBM,TMIX,SMIX,TN,SN,TERP)
*
*-----------------------------------------------------------------------
*
*Purpose:
* Compute the Lagrange interpolation factors (TERP) for temperature and
* dilution interpolation of cross sections. TRANSX 2.0 algorithm.
*
*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
* NGRO number of energy groups.
* NSUBM number of submaterials (number of temperature/dilution
* collocation points).
* TMIX temperature of each submaterial in the library.
* SMIX dilution of each submaterial in the library.
* The submaterials are ordered by decreasing dilution and
* then by increasing temperature.
* TN temperature of the isotope.
* SN dilution cross section in each energy group of the isotope.
* A value of 1.0E10 is used for infinite dilution.
*
*Parameters: output
* TERP Lagrange interpolation factors.
*
*-----------------------------------------------------------------------
*
*----
* SUBROUTINE ARGUMENTS
*----
INTEGER NGRO,NSUBM
REAL TMIX(NSUBM),SMIX(NSUBM),TN,SN(NGRO),TERP(NGRO,NSUBM)
REAL, ALLOCATABLE, DIMENSION(:) :: WORK
*----
* LOCAL VARIABLES
*----
PARAMETER (IOUT=6)
DOUBLE PRECISION CHECK
*
IF(NSUBM.EQ.1) THEN
DO 10 I=1,NGRO
TERP(I,1)=1.0
10 CONTINUE
RETURN
ENDIF
*----
* SCRATCH STORAGE ALLOCATION
*----
ALLOCATE(WORK(NSUBM))
*
DO 115 ISUBM=1,NSUBM
TMAT=TMIX(ISUBM)
SMAT=SMIX(ISUBM)
*----
* FIND TEMPERATURES VALUES
*----
NTEMP=0
DO 20 JSUBM=1,NSUBM
IF(SMIX(JSUBM).LT.0.9E10) GO TO 20
NTEMP=NTEMP+1
WORK(NTEMP)=TMIX(JSUBM)
20 CONTINUE
TERPT=0.0
*----
* LIMIT TEMPERATURE INTERPOLATION TO LINEAR IS TN A GRID TEMPERATURE?
*----
DO 30 ITMP=1,NTEMP
TT=WORK(ITMP)
IF(ABS(TN-TT).LT.1.E-3*TN+1.E-3) THEN
IF(ABS(TN-TMAT).LT.1.E-3*TMAT+1.E-3) TERPT=1.0
GO TO 70
ENDIF
30 CONTINUE
*----
* IF TEMP OUT OF RANGE USE ENDPOINTS
*----
IF((TN.GT.WORK(NTEMP).AND.ABS(TMAT-WORK(NTEMP)).LT.1.E-3*TMAT
> +1.E-3).OR.(TN.LT.WORK(1).AND.ABS(TMAT-WORK(1)).LT.
> 1.E-3*TMAT+1.E-3)) THEN
TERPT=1.0
GO TO 70
ENDIF
*----
* FIND BRACKETING TEMPS
*----
IF(NTEMP.EQ.1) THEN
TERPT=1.0
GO TO 70
ENDIF
DO 40 ITMP=1,NTEMP-1
IF(WORK(ITMP).LT.TN.AND.WORK(ITMP+1).GT.TN) THEN
ILOW=ITMP
IHIGH=ITMP+1
IF(ABS(TMAT-WORK(ITMP)).LT.1.E-3*TMAT+1.E-3.OR.ABS(TMAT-
> WORK(ITMP+1)).LT.1.E-3*TMAT+1.E-3) GO TO 50
ENDIF
40 CONTINUE
GO TO 70
50 TERPT=1.0
DO 60 ITMP=ILOW,IHIGH
TT=WORK(ITMP)
IF(ABS(TMAT-TT).LT.1.E-3*TMAT+1.E-3) GO TO 60
TERPT=TERPT*(TN-TT)/(TMAT-TT)
IF(ABS(TERPT).LT.1.E-3) GO TO 70
60 CONTINUE
*
*----
* FIND SIGMA-ZERO VALUES
*----
70 NTEMP=0
NSIGZ=0
DO 80 JSUBM=1,NSUBM
IF(SMIX(JSUBM).GE.0.9E10) NTEMP=NTEMP+1
IF(NTEMP.GT.1) GO TO 80
NSIGZ=NSIGZ+1
WORK(NSIGZ)=SMIX(JSUBM)
80 CONTINUE
*----
* FIND TERP FACTOR FOR SIGMA-ZERO
*----
DO 110 JJ=1,NGRO
TERPS=0.0
IF((SN(JJ).GE.WORK(1)).OR.(NSIGZ.EQ.1)) THEN
IF(SMAT.EQ.WORK(1)) TERPS=1.0
ELSE IF(SN(JJ).LE.WORK(NSIGZ)) THEN
IF(SMAT.EQ.WORK(NSIGZ)) TERPS=1.0
ELSE IF((SN(JJ).GT.WORK(2)).OR.(NSIGZ.EQ.2)) THEN
IF(SMAT.EQ.WORK(1)) TERPS=1.0-WORK(2)/SN(JJ)
IF(SMAT.EQ.WORK(2)) TERPS=WORK(2)/SN(JJ)
ELSE
DO 90 I=2,NSIGZ-1
IF(SN(JJ).LT.WORK(I+1)) GO TO 90
IF(SMAT.EQ.WORK(I+1)) TERPS=(ALOG10(WORK(I))-ALOG10(SN(JJ)))
X /(ALOG10(WORK(I))-ALOG10(WORK(I+1)))
IF(SMAT.EQ.WORK(I)) TERPS=(ALOG10(SN(JJ))-ALOG10(WORK(I+1)))
X /(ALOG10(WORK(I))-ALOG10(WORK(I+1)))
GO TO 100
90 CONTINUE
ENDIF
100 TERP(JJ,ISUBM)=TERPT*TERPS
110 CONTINUE
115 CONTINUE
*----
* CHECK FOR CONSISTENCY OF THE TERP FACTORS.
*----
DO 130 JJ=1,NGRO
CHECK=0.0D0
DO 120 ISUBM=1,NSUBM
CHECK=CHECK+TERP(JJ,ISUBM)
120 CONTINUE
IF(ABS(CHECK-1.0D0).GT.5.0D-3) THEN
WRITE (IOUT,200) JJ,CHECK,(TERP(JJ,ISUBM),ISUBM=1,NSUBM)
CALL XABORT('LIBTE2: INTERPOLATION FAILURE.')
ENDIF
130 CONTINUE
*----
* SCRATCH STORAGE DEALLOCATION
*----
DEALLOCATE(WORK)
RETURN
*
200 FORMAT (/51H LIBTE2: INCONSISTENT LAGRANGE INTERPOLATION FACTOR,
1 9H IN GROUP,I4,8H. CHECK=,1P,E13.3/9H FACTORS=,1P,9E13.3/
2 (9X,1P,9E13.3))
END
|