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 G4double momentum_mag_square = y[3]*y[3] + y[4]*y[4] + y[5]*y[5];
108 G4double inv_momentum_magnitude = 1.0 / std::sqrt( momentum_mag_square );
109
112
116
117 dydx[0] = y[3]*inv_momentum_magnitude;
118 dydx[1] = y[4]*inv_momentum_magnitude;
119 dydx[2] = y[5]*inv_momentum_magnitude;
120
121 dydx[3] = 0.;
122 dydx[4] = 0.;
123 dydx[5] = 0.;
124
125 G4double field[18] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
126
127 field[0] = Field[0];
128 field[1] = Field[1];
129 field[2] = Field[2];
130
131
132
134 {
136 {
137 dydx[3] += cof1*(y[4]*field[2] - y[5]*field[1]);
138 dydx[4] += cof1*(y[5]*field[0] - y[3]*field[2]);
139 dydx[5] += cof1*(y[3]*field[1] - y[4]*field[0]);
140 }
141 }
142
143
144
146 {
147 field[3] = Field[0];
148 field[4] = Field[1];
149 field[5] = Field[2];
150 }
151 else
152 {
153 field[3] = Field[3];
154 field[4] = Field[4];
155 field[5] = Field[5];
156 }
157
159 {
161 {
162 dydx[3] += cof1*cof2*field[3];
163 dydx[4] += cof1*cof2*field[4];
164 dydx[5] += cof1*cof2*field[5];
165 }
166 }
167
168
169
171 {
172 field[6] = Field[0];
173 field[7] = Field[1];
174 field[8] = Field[2];
175 }
176 else
177 {
178 field[6] = Field[6];
179 field[7] = Field[7];
180 field[8] = Field[8];
181 }
182
184 {
186 {
187 dydx[3] += field[6]*cof2*cof3/
c_light;
188 dydx[4] += field[7]*cof2*cof3/
c_light;
189 dydx[5] += field[8]*cof2*cof3/
c_light;
190 }
191 }
192
193
194
196 {
197 field[9] = Field[0];
198 field[10] = Field[1];
199 field[11] = Field[2];
200 field[12] = Field[3];
201 field[13] = Field[4];
202 field[14] = Field[5];
203 field[15] = Field[6];
204 field[16] = Field[7];
205 field[17] = Field[8];
206 }
207 else
208 {
209 field[9] = Field[9];
210 field[10] = Field[10];
211 field[11] = Field[11];
212 field[12] = Field[12];
213 field[13] = Field[13];
214 field[14] = Field[14];
215 field[15] = Field[15];
216 field[16] = Field[16];
217 field[17] = Field[17];
218 }
219
221 {
223 {
224 dydx[3] +=
magMoment*(y[9]*field[ 9]+y[10]*field[10]+y[11]*field[11])
225 *inv_momentum_magnitude*Energy;
226 dydx[4] +=
magMoment*(y[9]*field[12]+y[10]*field[13]+y[11]*field[14])
227 *inv_momentum_magnitude*Energy;
228 dydx[5] +=
magMoment*(y[9]*field[15]+y[10]*field[16]+y[11]*field[17])
229 *inv_momentum_magnitude*Energy;
230 }
231 }
232
233 dydx[6] = 0.;
234
235
236
237 dydx[7] = inverse_velocity;
238
240 {
241 dydx[ 8] = 0.;
242
243 dydx[ 9] = 0.;
244 dydx[10] = 0.;
245 dydx[11] = 0.;
246 }
247
249 {
252 {
254 BField = F;
255 }
256
259 {
261 EField = F;
262 }
263
265
267 u *= inv_momentum_magnitude;
268
272
274
276 if (
charge == 0.) pcharge = 1.;
278
280 if (
Spin.mag2() != 0.)
281 {
283 {
284 dSpin =
286 }
288 {
289 dSpin -= pcharge*
omegac*( uce*(u*(
Spin*EField) - EField*(Spin*u)) );
290
291
292
293 }
294 }
295
296 dydx[ 9] = dSpin.x();
297 dydx[10] = dSpin.y();
298 dydx[11] = dSpin.z();
299 }
300
301 return;
302}
