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
|
/*******************************************************************************
*
* McStas, neutron ray-tracing package
* Copyright 1997-2002, All rights reserved
* Risoe National Laboratory, Roskilde, Denmark
* Institut Laue Langevin, Grenoble, France
*
* Component: Source_Maxwell_3
*
* %I
* Written by: Kim Lefmann
* Date: March 2001
* Origin: Risoe
*
* Source with up to three Maxwellian distributions
*
* %D
* A parametrised continuous source for modelling a (cubic) source
* with (up to) 3 Maxwellian distributions.
* The source produces a continuous spectrum.
* The sampling of the neutrons is uniform in wavelength.
*
* Units of flux: neutrons/cm^2/second/ster
* (McStas units are in general neutrons/second)
*
* Example: PSI cold source T1=150.42 K / 2.51 AA I1 = 3.67 E11
* T2=38.74 K / 4.95 AA I2 = 3.64 E11
* T3=14.84 K / 9.5 AA I3 = 0.95 E11
*
* %P
* Input parameters:
*
* yheight: [m] Height of rectangular source
* xwidth: [m] Width of rectangular source
* Lmin: [AA] Lower edge of lambda distribution
* Lmax: [AA] Upper edge of lambda distribution
* lambda0: [AA] Mean wavelength of neutrons.
* dlambda: [AA] Wavelength spread of neutrons.
* target_index: [1] relative index of component to focus at, e.g. next is +1 this is used to compute 'dist' automatically.
* focus_xw: [m] Width of focusing rectangle
* focus_yh: [m] Height of focusing rectangle
* T1: [K] 1st temperature of thermal distribution
* I1: [1/(cm**2*st)] flux, 1 (in flux units, see above)
*
* Optional parameters:
* dist: [m] Distance from source to focusing rectangle; at (0,0,dist)
* T2: [K] 2nd temperature of thermal distribution
* T3: [K] 3nd temperature of - - -
* I2: [1/(cm**2*st)] flux, 2 (in flux units, see above)
* I3: [1/(cm**2*st)] flux, 3 - - -
* size: [m] Edge of cube shaped source (for backward compatibility)
*
* %E
*******************************************************************************/
DEFINE COMPONENT Source_Maxwell_3
SETTING PARAMETERS (size=0, yheight=0, xwidth=0, Lmin, Lmax, dist,
focus_xw, focus_yh,
T1, T2=300, T3=300, I1, I2=0, I3=0,
int target_index=+1,lambda0=0, dlambda=0)
/* Neutron parameters: (x,y,z,vx,vy,vz,t,sx,sy,sz,p) */
SHARE
%{
/* A normalised Maxwellian distribution : Integral over all l = 1 */
#pragma acc routine seq
double SM3_Maxwell(double l, double temp)
{
double a=949.0/temp;
return 2*a*a*exp(-a/(l*l))/(l*l*l*l*l);
}
%}
DECLARE
%{
double l_range;
double w_mult;
double w_source;
double h_source;
%}
INITIALIZE
%{
if (target_index && !dist)
{
Coords ToTarget;
double tx,ty,tz;
ToTarget = coords_sub(POS_A_COMP_INDEX(INDEX_CURRENT_COMP+target_index),POS_A_CURRENT_COMP);
ToTarget = rot_apply(ROT_A_CURRENT_COMP, ToTarget);
coords_get(ToTarget, &tx, &ty, &tz);
dist=sqrt(tx*tx+ty*ty+tz*tz);
}
if (size>0) {
w_source = h_source = size;
} else {
w_source = xwidth;
h_source = yheight;
}
if (lambda0) {
Lmin=lambda0-dlambda;
Lmax=lambda0+dlambda;
}
l_range = Lmax-Lmin;
w_mult = w_source*h_source*1.0e4; /* source area correction */
w_mult *= l_range; /* wavelength range correction */
w_mult *= 1.0/mcget_ncount(); /* correct for # neutron rays */
if (w_source <0 || h_source < 0 || Lmin <= 0 || Lmax <= 0 || dist <= 0 || T1 <= 0 || T2 <= 0|| T3 <= 0 || Lmax<=Lmin) {
printf("Source_Maxwell_3: %s: Error in input parameter values!\n"
"ERROR Exiting\n",
NAME_CURRENT_COMP);
exit(0);
}
%}
TRACE
%{
double v,tau_l,E,lambda,k,r,xf,yf,dx,dy,w_focus;
t=0;
z=0;
x = 0.5*w_source*randpm1();
y = 0.5*h_source*randpm1(); /* Choose initial position */
randvec_target_rect_real(&xf, &yf, &r, &w_focus,
0, 0, dist, focus_xw, focus_yh, ROT_A_CURRENT_COMP, x, y, z, 2);
dx = xf-x;
dy = yf-y;
r = sqrt(dx*dx+dy*dy+dist*dist);
lambda = Lmin+l_range*rand01(); /* Choose from uniform distribution */
k = 2*PI/lambda;
v = K2V*k;
vz = v*dist/r;
vy = v*dy/r;
vx = v*dx/r;
/* printf("pos0 (%g %g %g), pos1 (%g %g %g), r: %g, v (%g %g %g), v %g\n",
x,y,z,xf,yf,dist,r,vx,vy,vz, v);
printf("l %g, w_focus %g \n", lambda, w_focus); */
p *= w_mult*w_focus; /* Correct for target focusing etc */
p *= I1*SM3_Maxwell(lambda,T1)+I2*SM3_Maxwell(lambda,T2)+I3*SM3_Maxwell(lambda,T3);
/* Calculate true intensity */
%}
MCDISPLAY
%{
multiline(5, -(double)focus_xw/2.0, -(double)focus_yh/2.0, 0.0,
(double)focus_xw/2.0, -(double)focus_yh/2.0, 0.0,
(double)focus_xw/2.0, (double)focus_yh/2.0, 0.0,
-(double)focus_xw/2.0, (double)focus_yh/2.0, 0.0,
-(double)focus_xw/2.0, -(double)focus_yh/2.0, 0.0);
if (dist) {
dashed_line(0,0,0, -focus_xw/2,-focus_yh/2,dist, 4);
dashed_line(0,0,0, focus_xw/2,-focus_yh/2,dist, 4);
dashed_line(0,0,0, focus_xw/2, focus_yh/2,dist, 4);
dashed_line(0,0,0, -focus_xw/2, focus_yh/2,dist, 4);
}
%}
END
|
