]>
Commit | Line | Data |
---|---|---|
93d28366 | 1 | *$ CREATE FLUSCW.FOR |
2 | *COPY FLUSCW | |
3 | * * | |
4 | *=== fluscw ===========================================================* | |
5 | * * | |
6 | DOUBLE PRECISION FUNCTION FLUSCW | |
7 | #(IJ,PLA,TXX,TYY,TZZ,WEE,XX,YY,ZZ,NREG,IOLREG,LLO,NSURF) | |
8 | ||
9 | INCLUDE '(DBLPRC)' | |
10 | INCLUDE '(DIMPAR)' | |
11 | INCLUDE '(IOUNIT)' | |
12 | SAVE | |
13 | ||
14 | INCLUDE '(PAPROP)' | |
15 | INCLUDE '(USRBIN)' | |
16 | INCLUDE '(USRBDX)' | |
17 | INCLUDE '(USRTRC)' | |
18 | INCLUDE '(SCOHLP)' | |
19 | * | |
20 | *----------------------------------------------------------------------* | |
21 | * * | |
22 | * This functions returns neutron, proton and pion displacement damage * | |
23 | * weight factors. * | |
24 | * * | |
25 | *----------------------------------------------------------------------* | |
26 | * * | |
27 | * Input variables: * | |
28 | * * | |
29 | * Ij = (generalized) particle code * | |
30 | * Pla = particle momentum (if > 0), or kinetic energy (if <0 )* | |
31 | * Txx,yy,zz = particle direction cosines * | |
32 | * Wee = particle weight * | |
33 | * Xx,Yy,Zz = position * | |
34 | * Nreg = (new) region number * | |
35 | * Iolreg = (old) region number * | |
36 | * Llo = particle generation * | |
37 | * Nsurf = transport flag (ignore!) * | |
38 | * * | |
39 | * Output variables: * | |
40 | * * | |
41 | * Fluscw = factor the scored amount will be multiplied by * | |
42 | * Lsczer = logical flag, if true no amount will be scored * | |
43 | * regardless of Fluscw * | |
44 | * * | |
45 | * Useful variables (common SCOHLP): * | |
46 | * * | |
47 | * Flux like binnings/estimators (Fluscw): * | |
48 | * ISCRNG = 1 --> Boundary crossing estimator * | |
49 | * ISCRNG = 2 --> Track length binning * | |
50 | * ISCRNG = 3 --> Track length estimator * | |
51 | * ISCRNG = 4 --> Collision density estimator * | |
52 | * ISCRNG = 5 --> Yield estimator * | |
53 | * JSCRNG = # of the binning/estimator * | |
54 | * * | |
55 | *----------------------------------------------------------------------* | |
56 | * | |
57 | LOGICAL LFIRST | |
58 | ||
59 | ||
60 | DATA LFIRST /.TRUE./ | |
61 | * | |
62 | * the default is unit weight | |
63 | FLUSCW = ONEONE | |
64 | LSCZER = .FALSE. | |
65 | * | |
66 | * skip all scorings other than tracklength | |
67 | IF (ISCRNG.NE.2) RETURN | |
68 | * skip all particles other than protons, pions and neutrons | |
69 | IF ((IJ.NE.1).AND.(IJ.NE.13).AND.(IJ.NE.14).AND.(IJ.NE.8)) RETURN | |
70 | * | |
71 | * read displacement damage factors | |
72 | IF (LFIRST) THEN | |
73 | WRITE(LUNOUT,1000) | |
74 | 1000 FORMAT(1X,'FLUSCW: direct conversion of neutron fluence to', | |
75 | & ' 1 MeV neutron equivalent requested') | |
76 | WRITE(LUNOUT,*) ' neutrons' | |
77 | OPEN(19,FILE='disdam_n.dat',STATUS='UNKNOWN') | |
78 | CALL SKIP(19,7) | |
79 | CALL RDXSC(19,LUNOUT,IDNDAM,-3) | |
80 | CLOSE(19) | |
81 | WRITE(LUNOUT,*) ' protons ' | |
82 | OPEN(19,FILE='disdam_p.dat',STATUS='UNKNOWN') | |
83 | CALL SKIP(19,7) | |
84 | CALL RDXSC(19,LUNOUT,IDPDAM,-2) | |
85 | CLOSE(19) | |
86 | WRITE(LUNOUT,*) ' pions ' | |
87 | OPEN(19,FILE='disdam_pi.dat',STATUS='UNKNOWN') | |
88 | CALL SKIP(19,7) | |
89 | CALL RDXSC(19,LUNOUT,IDODAM,-2) | |
90 | CLOSE(19) | |
91 | LFIRST = .FALSE. | |
92 | ENDIF | |
93 | * | |
94 | * should be always called with ekin ( pla < 0 ), | |
95 | * but we leave the check for the moment.. | |
96 | IF (PLA.LT.0.0D0) THEN | |
97 | EKIN = ABS(PLA) | |
98 | ELSEIF (PLA.GT.0.0D0) THEN | |
99 | EKIN = SQRT(PLA**2+AM(IJ)**2)-AM(IJ) | |
100 | ELSE | |
101 | RETURN | |
102 | ENDIF | |
103 | * | |
104 | * calculate the weight | |
105 | IF (IJ.EQ.1) THEN | |
106 | FLUSCW = XSECT(IDPDAM,EKIN,0) | |
107 | ELSEIF ((IJ.EQ.13).OR.(IJ.EQ.14)) THEN | |
108 | FLUSCW = XSECT(IDODAM,EKIN,0) | |
109 | ELSEIF (IJ.EQ.8) THEN | |
110 | FLUSCW = XSECT(IDNDAM,EKIN,0) | |
111 | ELSE | |
112 | STOP ' FLUSCW: inconsistent IJ ' | |
113 | ENDIF | |
114 | ||
115 | RETURN | |
116 | END | |
117 | * | |
118 | *===rdxsc==============================================================* | |
119 | * | |
120 | SUBROUTINE RDXSC(LINP,LCHK,IDXSEC,MODE) | |
121 | ||
122 | ************************************************************************ | |
123 | * LINP logical input unit * | |
124 | * LCHK > 0 cross section data are plotted to unit LCHK * | |
125 | * in equidistant log. binning using double-log. * | |
126 | * interpolation * | |
127 | * (otherwise they are not plotted) * | |
128 | * IDXSEC index of stored cross section in common * | |
129 | * |MODE| = 1 cross section data given as histogram * | |
130 | * i.e. E_lo(i) sigma(i) * | |
131 | * E_hi(i) sigma(i) * | |
132 | * with increasing energy * | |
133 | * = 2 cross section data given for bin averages with * | |
134 | * increasing energy * | |
135 | * = 3 cross section data given for bin averages with * | |
136 | * decreasing energy * | |
137 | * if MODE < 0 the energy unit is assumed to be MeV (otherwise GeV) * | |
138 | ************************************************************************ | |
139 | ||
140 | IMPLICIT DOUBLE PRECISION (A-H,O-Z) | |
141 | SAVE | |
142 | ||
143 | PARAMETER (MAXSDT = 10, | |
144 | & MAXSBI = 1400) | |
145 | COMMON /CSECT/ XSEAV(MAXSDT,MAXSBI),XS(MAXSDT,MAXSBI), | |
146 | & NXSBIN(MAXSDT),NXSDT | |
147 | DIMENSION TMPXS(2,MAXSBI) | |
148 | ||
149 | LOGICAL LFIRST | |
150 | DATA LFIRST /.TRUE./ | |
151 | ||
152 | IF (LFIRST) THEN | |
153 | NXSDT = 0 | |
154 | DO 1 IDT=1,MAXSDT | |
155 | NXSBIN(IDT) = 0 | |
156 | DO 2 IBIN=1,MAXSBI | |
157 | XSEAV(IDT,IBIN) = 0.0D0 | |
158 | XS(IDT,IBIN) = 0.0D0 | |
159 | 2 CONTINUE | |
160 | 1 CONTINUE | |
161 | LFIRST = .FALSE. | |
162 | ENDIF | |
163 | ||
164 | NXSDT = NXSDT+1 | |
165 | ||
166 | NEBIN = 0 | |
167 | 10 CONTINUE | |
168 | NEBIN = NEBIN+1 | |
169 | IF (NEBIN.GT.MAXSBI) STOP ' nebin > maxsbi !!!' | |
170 | IF (IABS(MODE).EQ.1) THEN | |
171 | READ(LINP,*,END=9000) ELOW,XS(NXSDT,NEBIN) | |
172 | READ(LINP,*) ELHI,XS(NXSDT,NEBIN) | |
173 | XSEAV(NXSDT,NEBIN) = SQRT(ELOW*ELHI) | |
174 | IF (MODE.LT.0) XSEAV(NXSDT,NEBIN) = 1.D-3*XSEAV(NXSDT,NEBIN) | |
175 | ELSEIF (IABS(MODE).EQ.2) THEN | |
176 | READ(LINP,*,END=9000) XSEAV(NXSDT,NEBIN),XS(NXSDT,NEBIN) | |
177 | IF (MODE.LT.0) XSEAV(NXSDT,NEBIN) = 1.D-3*XSEAV(NXSDT,NEBIN) | |
178 | ELSEIF (IABS(MODE).EQ.3) THEN | |
179 | READ(LINP,*,END=9000) TMPXS(1,NEBIN),TMPXS(2,NEBIN) | |
180 | IF (MODE.LT.0) TMPXS(1,NEBIN) = 1.D-3*TMPXS(1,NEBIN) | |
181 | ELSE | |
182 | STOP ' RDXSC: unsupported mode ! ' | |
183 | ENDIF | |
184 | GOTO 10 | |
185 | 9000 CONTINUE | |
186 | ||
187 | IF (IABS(MODE).EQ.3) THEN | |
188 | DO 3 I=1,NEBIN-1 | |
189 | IDX = NEBIN-I | |
190 | XSEAV(NXSDT,IDX) = TMPXS(1,I) | |
191 | XS(NXSDT,IDX) = TMPXS(2,I) | |
192 | 3 CONTINUE | |
193 | ENDIF | |
194 | ||
195 | NXSBIN(NXSDT) = NEBIN-1 | |
196 | IDXSEC = NXSDT | |
197 | ||
198 | IF (LCHK.GT.0) THEN | |
199 | AELO = LOG10(XSEAV(IDXSEC,1)) | |
200 | AEHI = LOG10(XSEAV(IDXSEC,NXSBIN(IDXSEC))) | |
201 | DAE = (AEHI-AELO)/DBLE(NXSBIN(IDXSEC)) | |
202 | DO 4 I=1,NXSBIN(IDXSEC)+1 | |
203 | AE = AELO+DBLE(I-1)*DAE | |
204 | E = 10.0D0**AE | |
205 | WRITE(LCHK,'(1X,2E15.5)') E,XSECT(IDXSEC,E,0) | |
206 | 4 CONTINUE | |
207 | ENDIF | |
208 | ||
209 | RETURN | |
210 | END | |
211 | * | |
212 | *===xsect==============================================================* | |
213 | * | |
214 | DOUBLE PRECISION FUNCTION XSECT(IDXSEC,E,MODE) | |
215 | ||
216 | ************************************************************************ | |
217 | * IDXSEC index of stored cross section in common * | |
218 | * E energy * | |
219 | ************************************************************************ | |
220 | ||
221 | IMPLICIT DOUBLE PRECISION (A-H,O-Z) | |
222 | SAVE | |
223 | ||
224 | PARAMETER (MAXSDT = 10, | |
225 | & MAXSBI = 1400) | |
226 | COMMON /CSECT/ XSEAV(MAXSDT,MAXSBI),XS(MAXSDT,MAXSBI), | |
227 | & NXSBIN(MAXSDT),NXSDT | |
228 | ||
229 | IF (E.LE.XSEAV(IDXSEC,1)) THEN | |
230 | CS = 0.0D0 | |
231 | ELSEIF (E.GT.XSEAV(IDXSEC,NXSBIN(IDXSEC))) THEN | |
232 | N = NXSBIN(IDXSEC)-1 | |
233 | CS = XSCINT(IDXSEC,E,N) | |
234 | ELSE | |
235 | DO 1 J=1,NXSBIN(IDXSEC)-1 | |
236 | IF ((E.GT.XSEAV(IDXSEC,J)).AND.(E.LE.XSEAV(IDXSEC,J+1))) | |
237 | & THEN | |
238 | CS = XSCINT(IDXSEC,E,J) | |
239 | GOTO 2 | |
240 | ENDIF | |
241 | 1 CONTINUE | |
242 | STOP ' xsection value not found ' | |
243 | 2 CONTINUE | |
244 | ENDIF | |
245 | XSECT = CS | |
246 | ||
247 | RETURN | |
248 | END | |
249 | * | |
250 | *===xscint=============================================================* | |
251 | * | |
252 | DOUBLE PRECISION FUNCTION XSCINT(IDXSEC,E,J) | |
253 | ||
254 | IMPLICIT DOUBLE PRECISION (A-H,O-Z) | |
255 | SAVE | |
256 | ||
257 | PARAMETER (MAXSDT = 10, | |
258 | & MAXSBI = 1400) | |
259 | COMMON /CSECT/ XSEAV(MAXSDT,MAXSBI),XS(MAXSDT,MAXSBI), | |
260 | & NXSBIN(MAXSDT),NXSDT | |
261 | ||
262 | FAC = (LOG10(E) -LOG10(XSEAV(IDXSEC,J)))/ | |
263 | & (LOG10(XSEAV(IDXSEC,J+1))-LOG10(XSEAV(IDXSEC,J))) | |
264 | XSCINT = LOG10(XS(IDXSEC,J)) | |
265 | & +(LOG10(XS(IDXSEC,J+1))-LOG10(XS(IDXSEC,J)))*FAC | |
266 | XSCINT = 10**(XSCINT) | |
267 | ||
268 | RETURN | |
269 | END | |
270 | * | |
271 | *===skip===============================================================* | |
272 | * | |
273 | SUBROUTINE SKIP(LUNIT,NSKIP) | |
274 | ||
275 | IMPLICIT DOUBLE PRECISION (A-H,O-Z) | |
276 | SAVE | |
277 | PARAMETER (LOUT=6,LLOOK=9) | |
278 | ||
279 | CHARACTER*132 ACARD | |
280 | ||
281 | IF (NSKIP.EQ.0) RETURN | |
282 | ||
283 | DO 1 K=1,NSKIP | |
284 | READ(LUNIT,'(A132)') ACARD | |
285 | 1 CONTINUE | |
286 | ||
287 | RETURN | |
288 | END |