CLHEP 2.4.6.4
C++ Class Library for High Energy Physics
Loading...
Searching...
No Matches
MatrixInvert.cc
Go to the documentation of this file.
1// -*- C++ -*-
2//
3// This file is a part of the CLHEP - a Class Library for High Energy Physics.
4//
5// This software written by Mark Fischler and Steven Haywood
6//
7
8#include <string.h>
9
10#include "CLHEP/Matrix/defs.h"
11#include "CLHEP/Matrix/Matrix.h"
12
13namespace CLHEP {
14
15// Aij are indices for a 6x6 matrix.
16// Mij are indices for a 5x5 matrix.
17// Fij are indices for a 4x4 matrix.
18
19#define A00 0
20#define A01 1
21#define A02 2
22#define A03 3
23#define A04 4
24#define A05 5
25
26#define A10 6
27#define A11 7
28#define A12 8
29#define A13 9
30#define A14 10
31#define A15 11
32
33#define A20 12
34#define A21 13
35#define A22 14
36#define A23 15
37#define A24 16
38#define A25 17
39
40#define A30 18
41#define A31 19
42#define A32 20
43#define A33 21
44#define A34 22
45#define A35 23
46
47#define A40 24
48#define A41 25
49#define A42 26
50#define A43 27
51#define A44 28
52#define A45 29
53
54#define A50 30
55#define A51 31
56#define A52 32
57#define A53 33
58#define A54 34
59#define A55 35
60
61#define M00 0
62#define M01 1
63#define M02 2
64#define M03 3
65#define M04 4
66
67#define M10 5
68#define M11 6
69#define M12 7
70#define M13 8
71#define M14 9
72
73#define M20 10
74#define M21 11
75#define M22 12
76#define M23 13
77#define M24 14
78
79#define M30 15
80#define M31 16
81#define M32 17
82#define M33 18
83#define M34 19
84
85#define M40 20
86#define M41 21
87#define M42 22
88#define M43 23
89#define M44 24
90
91#define F00 0
92#define F01 1
93#define F02 2
94#define F03 3
95
96#define F10 4
97#define F11 5
98#define F12 6
99#define F13 7
100
101#define F20 8
102#define F21 9
103#define F22 10
104#define F23 11
105
106#define F30 12
107#define F31 13
108#define F32 14
109#define F33 15
110
111␌
112void HepMatrix::invertHaywood4 (int & ifail) {
113
114 ifail = 0;
115
116 // Find all NECESSARY 2x2 dets: (18 of them)
117
118 double Det2_12_01 = m[F10]*m[F21] - m[F11]*m[F20];
119 double Det2_12_02 = m[F10]*m[F22] - m[F12]*m[F20];
120 double Det2_12_03 = m[F10]*m[F23] - m[F13]*m[F20]; //
121 double Det2_12_13 = m[F11]*m[F23] - m[F13]*m[F21]; //
122 double Det2_12_23 = m[F12]*m[F23] - m[F13]*m[F22]; //
123 double Det2_12_12 = m[F11]*m[F22] - m[F12]*m[F21];
124 double Det2_13_01 = m[F10]*m[F31] - m[F11]*m[F30];
125 double Det2_13_02 = m[F10]*m[F32] - m[F12]*m[F30];
126 double Det2_13_03 = m[F10]*m[F33] - m[F13]*m[F30];
127 double Det2_13_12 = m[F11]*m[F32] - m[F12]*m[F31];
128 double Det2_13_13 = m[F11]*m[F33] - m[F13]*m[F31];
129 double Det2_13_23 = m[F12]*m[F33] - m[F13]*m[F32]; //
130 double Det2_23_01 = m[F20]*m[F31] - m[F21]*m[F30];
131 double Det2_23_02 = m[F20]*m[F32] - m[F22]*m[F30];
132 double Det2_23_03 = m[F20]*m[F33] - m[F23]*m[F30];
133 double Det2_23_12 = m[F21]*m[F32] - m[F22]*m[F31];
134 double Det2_23_13 = m[F21]*m[F33] - m[F23]*m[F31];
135 double Det2_23_23 = m[F22]*m[F33] - m[F23]*m[F32];
136
137 // Find all NECESSARY 3x3 dets: (16 of them)
138
139 double Det3_012_012 = m[F00]*Det2_12_12 - m[F01]*Det2_12_02
140 + m[F02]*Det2_12_01;
141 double Det3_012_013 = m[F00]*Det2_12_13 - m[F01]*Det2_12_03
142 + m[F03]*Det2_12_01; //
143 double Det3_012_023 = m[F00]*Det2_12_23 - m[F02]*Det2_12_03
144 + m[F03]*Det2_12_02; //
145 double Det3_012_123 = m[F01]*Det2_12_23 - m[F02]*Det2_12_13
146 + m[F03]*Det2_12_12; //
147 double Det3_013_012 = m[F00]*Det2_13_12 - m[F01]*Det2_13_02
148 + m[F02]*Det2_13_01;
149 double Det3_013_013 = m[F00]*Det2_13_13 - m[F01]*Det2_13_03
150 + m[F03]*Det2_13_01;
151 double Det3_013_023 = m[F00]*Det2_13_23 - m[F02]*Det2_13_03
152 + m[F03]*Det2_13_02; //
153 double Det3_013_123 = m[F01]*Det2_13_23 - m[F02]*Det2_13_13
154 + m[F03]*Det2_13_12; //
155 double Det3_023_012 = m[F00]*Det2_23_12 - m[F01]*Det2_23_02
156 + m[F02]*Det2_23_01;
157 double Det3_023_013 = m[F00]*Det2_23_13 - m[F01]*Det2_23_03
158 + m[F03]*Det2_23_01;
159 double Det3_023_023 = m[F00]*Det2_23_23 - m[F02]*Det2_23_03
160 + m[F03]*Det2_23_02;
161 double Det3_023_123 = m[F01]*Det2_23_23 - m[F02]*Det2_23_13
162 + m[F03]*Det2_23_12; //
163 double Det3_123_012 = m[F10]*Det2_23_12 - m[F11]*Det2_23_02
164 + m[F12]*Det2_23_01;
165 double Det3_123_013 = m[F10]*Det2_23_13 - m[F11]*Det2_23_03
166 + m[F13]*Det2_23_01;
167 double Det3_123_023 = m[F10]*Det2_23_23 - m[F12]*Det2_23_03
168 + m[F13]*Det2_23_02;
169 double Det3_123_123 = m[F11]*Det2_23_23 - m[F12]*Det2_23_13
170 + m[F13]*Det2_23_12;
171
172 // Find the 4x4 det:
173
174 double det = m[F00]*Det3_123_123
175 - m[F01]*Det3_123_023
176 + m[F02]*Det3_123_013
177 - m[F03]*Det3_123_012;
178
179 if ( det == 0 ) {
180#ifdef SINGULAR_DIAGNOSTICS
181 std::cerr << "Kramer's rule inversion of a singular 4x4 matrix: "
182 << *this << "\n";
183#endif
184 ifail = 1;
185 return;
186 }
187
188 double oneOverDet = 1.0/det;
189 double mn1OverDet = - oneOverDet;
190
191 m[F00] = Det3_123_123 * oneOverDet;
192 m[F01] = Det3_023_123 * mn1OverDet;
193 m[F02] = Det3_013_123 * oneOverDet;
194 m[F03] = Det3_012_123 * mn1OverDet;
195
196 m[F10] = Det3_123_023 * mn1OverDet;
197 m[F11] = Det3_023_023 * oneOverDet;
198 m[F12] = Det3_013_023 * mn1OverDet;
199 m[F13] = Det3_012_023 * oneOverDet;
200
201 m[F20] = Det3_123_013 * oneOverDet;
202 m[F21] = Det3_023_013 * mn1OverDet;
203 m[F22] = Det3_013_013 * oneOverDet;
204 m[F23] = Det3_012_013 * mn1OverDet;
205
206 m[F30] = Det3_123_012 * mn1OverDet;
207 m[F31] = Det3_023_012 * oneOverDet;
208 m[F32] = Det3_013_012 * mn1OverDet;
209 m[F33] = Det3_012_012 * oneOverDet;
210
211 return;
212}
213
214
215␌
216void HepMatrix::invertHaywood5 (int & ifail) {
217
218 ifail = 0;
219
220 // Find all NECESSARY 2x2 dets: (30 of them)
221
222 double Det2_23_01 = m[M20]*m[M31] - m[M21]*m[M30];
223 double Det2_23_02 = m[M20]*m[M32] - m[M22]*m[M30];
224 double Det2_23_03 = m[M20]*m[M33] - m[M23]*m[M30];
225 double Det2_23_04 = m[M20]*m[M34] - m[M24]*m[M30];
226 double Det2_23_12 = m[M21]*m[M32] - m[M22]*m[M31]; //
227 double Det2_23_13 = m[M21]*m[M33] - m[M23]*m[M31];
228 double Det2_23_14 = m[M21]*m[M34] - m[M24]*m[M31]; //
229 double Det2_23_23 = m[M22]*m[M33] - m[M23]*m[M32];
230 double Det2_23_24 = m[M22]*m[M34] - m[M24]*m[M32]; //
231 double Det2_23_34 = m[M23]*m[M34] - m[M24]*m[M33]; //
232 double Det2_24_01 = m[M20]*m[M41] - m[M21]*m[M40];
233 double Det2_24_02 = m[M20]*m[M42] - m[M22]*m[M40];
234 double Det2_24_03 = m[M20]*m[M43] - m[M23]*m[M40];
235 double Det2_24_04 = m[M20]*m[M44] - m[M24]*m[M40];
236 double Det2_24_12 = m[M21]*m[M42] - m[M22]*m[M41];
237 double Det2_24_13 = m[M21]*m[M43] - m[M23]*m[M41];
238 double Det2_24_14 = m[M21]*m[M44] - m[M24]*m[M41];
239 double Det2_24_23 = m[M22]*m[M43] - m[M23]*m[M42];
240 double Det2_24_24 = m[M22]*m[M44] - m[M24]*m[M42];
241 double Det2_24_34 = m[M23]*m[M44] - m[M24]*m[M43]; //
242 double Det2_34_01 = m[M30]*m[M41] - m[M31]*m[M40];
243 double Det2_34_02 = m[M30]*m[M42] - m[M32]*m[M40];
244 double Det2_34_03 = m[M30]*m[M43] - m[M33]*m[M40];
245 double Det2_34_04 = m[M30]*m[M44] - m[M34]*m[M40];
246 double Det2_34_12 = m[M31]*m[M42] - m[M32]*m[M41];
247 double Det2_34_13 = m[M31]*m[M43] - m[M33]*m[M41];
248 double Det2_34_14 = m[M31]*m[M44] - m[M34]*m[M41];
249 double Det2_34_23 = m[M32]*m[M43] - m[M33]*m[M42];
250 double Det2_34_24 = m[M32]*m[M44] - m[M34]*m[M42];
251 double Det2_34_34 = m[M33]*m[M44] - m[M34]*m[M43];
252
253 // Find all NECESSARY 3x3 dets: (40 of them)
254
255 double Det3_123_012 = m[M10]*Det2_23_12 - m[M11]*Det2_23_02
256 + m[M12]*Det2_23_01;
257 double Det3_123_013 = m[M10]*Det2_23_13 - m[M11]*Det2_23_03
258 + m[M13]*Det2_23_01;
259 double Det3_123_014 = m[M10]*Det2_23_14 - m[M11]*Det2_23_04
260 + m[M14]*Det2_23_01; //
261 double Det3_123_023 = m[M10]*Det2_23_23 - m[M12]*Det2_23_03
262 + m[M13]*Det2_23_02;
263 double Det3_123_024 = m[M10]*Det2_23_24 - m[M12]*Det2_23_04
264 + m[M14]*Det2_23_02; //
265 double Det3_123_034 = m[M10]*Det2_23_34 - m[M13]*Det2_23_04
266 + m[M14]*Det2_23_03; //
267 double Det3_123_123 = m[M11]*Det2_23_23 - m[M12]*Det2_23_13
268 + m[M13]*Det2_23_12;
269 double Det3_123_124 = m[M11]*Det2_23_24 - m[M12]*Det2_23_14
270 + m[M14]*Det2_23_12; //
271 double Det3_123_134 = m[M11]*Det2_23_34 - m[M13]*Det2_23_14
272 + m[M14]*Det2_23_13; //
273 double Det3_123_234 = m[M12]*Det2_23_34 - m[M13]*Det2_23_24
274 + m[M14]*Det2_23_23; //
275 double Det3_124_012 = m[M10]*Det2_24_12 - m[M11]*Det2_24_02
276 + m[M12]*Det2_24_01;
277 double Det3_124_013 = m[M10]*Det2_24_13 - m[M11]*Det2_24_03
278 + m[M13]*Det2_24_01;
279 double Det3_124_014 = m[M10]*Det2_24_14 - m[M11]*Det2_24_04
280 + m[M14]*Det2_24_01;
281 double Det3_124_023 = m[M10]*Det2_24_23 - m[M12]*Det2_24_03
282 + m[M13]*Det2_24_02;
283 double Det3_124_024 = m[M10]*Det2_24_24 - m[M12]*Det2_24_04
284 + m[M14]*Det2_24_02;
285 double Det3_124_034 = m[M10]*Det2_24_34 - m[M13]*Det2_24_04
286 + m[M14]*Det2_24_03; //
287 double Det3_124_123 = m[M11]*Det2_24_23 - m[M12]*Det2_24_13
288 + m[M13]*Det2_24_12;
289 double Det3_124_124 = m[M11]*Det2_24_24 - m[M12]*Det2_24_14
290 + m[M14]*Det2_24_12;
291 double Det3_124_134 = m[M11]*Det2_24_34 - m[M13]*Det2_24_14
292 + m[M14]*Det2_24_13; //
293 double Det3_124_234 = m[M12]*Det2_24_34 - m[M13]*Det2_24_24
294 + m[M14]*Det2_24_23; //
295 double Det3_134_012 = m[M10]*Det2_34_12 - m[M11]*Det2_34_02
296 + m[M12]*Det2_34_01;
297 double Det3_134_013 = m[M10]*Det2_34_13 - m[M11]*Det2_34_03
298 + m[M13]*Det2_34_01;
299 double Det3_134_014 = m[M10]*Det2_34_14 - m[M11]*Det2_34_04
300 + m[M14]*Det2_34_01;
301 double Det3_134_023 = m[M10]*Det2_34_23 - m[M12]*Det2_34_03
302 + m[M13]*Det2_34_02;
303 double Det3_134_024 = m[M10]*Det2_34_24 - m[M12]*Det2_34_04
304 + m[M14]*Det2_34_02;
305 double Det3_134_034 = m[M10]*Det2_34_34 - m[M13]*Det2_34_04
306 + m[M14]*Det2_34_03;
307 double Det3_134_123 = m[M11]*Det2_34_23 - m[M12]*Det2_34_13
308 + m[M13]*Det2_34_12;
309 double Det3_134_124 = m[M11]*Det2_34_24 - m[M12]*Det2_34_14
310 + m[M14]*Det2_34_12;
311 double Det3_134_134 = m[M11]*Det2_34_34 - m[M13]*Det2_34_14
312 + m[M14]*Det2_34_13;
313 double Det3_134_234 = m[M12]*Det2_34_34 - m[M13]*Det2_34_24
314 + m[M14]*Det2_34_23; //
315 double Det3_234_012 = m[M20]*Det2_34_12 - m[M21]*Det2_34_02
316 + m[M22]*Det2_34_01;
317 double Det3_234_013 = m[M20]*Det2_34_13 - m[M21]*Det2_34_03
318 + m[M23]*Det2_34_01;
319 double Det3_234_014 = m[M20]*Det2_34_14 - m[M21]*Det2_34_04
320 + m[M24]*Det2_34_01;
321 double Det3_234_023 = m[M20]*Det2_34_23 - m[M22]*Det2_34_03
322 + m[M23]*Det2_34_02;
323 double Det3_234_024 = m[M20]*Det2_34_24 - m[M22]*Det2_34_04
324 + m[M24]*Det2_34_02;
325 double Det3_234_034 = m[M20]*Det2_34_34 - m[M23]*Det2_34_04
326 + m[M24]*Det2_34_03;
327 double Det3_234_123 = m[M21]*Det2_34_23 - m[M22]*Det2_34_13
328 + m[M23]*Det2_34_12;
329 double Det3_234_124 = m[M21]*Det2_34_24 - m[M22]*Det2_34_14
330 + m[M24]*Det2_34_12;
331 double Det3_234_134 = m[M21]*Det2_34_34 - m[M23]*Det2_34_14
332 + m[M24]*Det2_34_13;
333 double Det3_234_234 = m[M22]*Det2_34_34 - m[M23]*Det2_34_24
334 + m[M24]*Det2_34_23;
335
336 // Find all NECESSARY 4x4 dets: (25 of them)
337
338 double Det4_0123_0123 = m[M00]*Det3_123_123 - m[M01]*Det3_123_023
339 + m[M02]*Det3_123_013 - m[M03]*Det3_123_012;
340 double Det4_0123_0124 = m[M00]*Det3_123_124 - m[M01]*Det3_123_024
341 + m[M02]*Det3_123_014 - m[M04]*Det3_123_012; //
342 double Det4_0123_0134 = m[M00]*Det3_123_134 - m[M01]*Det3_123_034
343 + m[M03]*Det3_123_014 - m[M04]*Det3_123_013; //
344 double Det4_0123_0234 = m[M00]*Det3_123_234 - m[M02]*Det3_123_034
345 + m[M03]*Det3_123_024 - m[M04]*Det3_123_023; //
346 double Det4_0123_1234 = m[M01]*Det3_123_234 - m[M02]*Det3_123_134
347 + m[M03]*Det3_123_124 - m[M04]*Det3_123_123; //
348 double Det4_0124_0123 = m[M00]*Det3_124_123 - m[M01]*Det3_124_023
349 + m[M02]*Det3_124_013 - m[M03]*Det3_124_012;
350 double Det4_0124_0124 = m[M00]*Det3_124_124 - m[M01]*Det3_124_024
351 + m[M02]*Det3_124_014 - m[M04]*Det3_124_012;
352 double Det4_0124_0134 = m[M00]*Det3_124_134 - m[M01]*Det3_124_034
353 + m[M03]*Det3_124_014 - m[M04]*Det3_124_013; //
354 double Det4_0124_0234 = m[M00]*Det3_124_234 - m[M02]*Det3_124_034
355 + m[M03]*Det3_124_024 - m[M04]*Det3_124_023; //
356 double Det4_0124_1234 = m[M01]*Det3_124_234 - m[M02]*Det3_124_134
357 + m[M03]*Det3_124_124 - m[M04]*Det3_124_123; //
358 double Det4_0134_0123 = m[M00]*Det3_134_123 - m[M01]*Det3_134_023
359 + m[M02]*Det3_134_013 - m[M03]*Det3_134_012;
360 double Det4_0134_0124 = m[M00]*Det3_134_124 - m[M01]*Det3_134_024
361 + m[M02]*Det3_134_014 - m[M04]*Det3_134_012;
362 double Det4_0134_0134 = m[M00]*Det3_134_134 - m[M01]*Det3_134_034
363 + m[M03]*Det3_134_014 - m[M04]*Det3_134_013;
364 double Det4_0134_0234 = m[M00]*Det3_134_234 - m[M02]*Det3_134_034
365 + m[M03]*Det3_134_024 - m[M04]*Det3_134_023; //
366 double Det4_0134_1234 = m[M01]*Det3_134_234 - m[M02]*Det3_134_134
367 + m[M03]*Det3_134_124 - m[M04]*Det3_134_123; //
368 double Det4_0234_0123 = m[M00]*Det3_234_123 - m[M01]*Det3_234_023
369 + m[M02]*Det3_234_013 - m[M03]*Det3_234_012;
370 double Det4_0234_0124 = m[M00]*Det3_234_124 - m[M01]*Det3_234_024
371 + m[M02]*Det3_234_014 - m[M04]*Det3_234_012;
372 double Det4_0234_0134 = m[M00]*Det3_234_134 - m[M01]*Det3_234_034
373 + m[M03]*Det3_234_014 - m[M04]*Det3_234_013;
374 double Det4_0234_0234 = m[M00]*Det3_234_234 - m[M02]*Det3_234_034
375 + m[M03]*Det3_234_024 - m[M04]*Det3_234_023;
376 double Det4_0234_1234 = m[M01]*Det3_234_234 - m[M02]*Det3_234_134
377 + m[M03]*Det3_234_124 - m[M04]*Det3_234_123; //
378 double Det4_1234_0123 = m[M10]*Det3_234_123 - m[M11]*Det3_234_023
379 + m[M12]*Det3_234_013 - m[M13]*Det3_234_012;
380 double Det4_1234_0124 = m[M10]*Det3_234_124 - m[M11]*Det3_234_024
381 + m[M12]*Det3_234_014 - m[M14]*Det3_234_012;
382 double Det4_1234_0134 = m[M10]*Det3_234_134 - m[M11]*Det3_234_034
383 + m[M13]*Det3_234_014 - m[M14]*Det3_234_013;
384 double Det4_1234_0234 = m[M10]*Det3_234_234 - m[M12]*Det3_234_034
385 + m[M13]*Det3_234_024 - m[M14]*Det3_234_023;
386 double Det4_1234_1234 = m[M11]*Det3_234_234 - m[M12]*Det3_234_134
387 + m[M13]*Det3_234_124 - m[M14]*Det3_234_123;
388
389 // Find the 5x5 det:
390
391 double det = m[M00]*Det4_1234_1234
392 - m[M01]*Det4_1234_0234
393 + m[M02]*Det4_1234_0134
394 - m[M03]*Det4_1234_0124
395 + m[M04]*Det4_1234_0123;
396
397 if ( det == 0 ) {
398#ifdef SINGULAR_DIAGNOSTICS
399 std::cerr << "Kramer's rule inversion of a singular 5x5 matrix: "
400 << *this << "\n";
401#endif
402 ifail = 1;
403 return;
404 }
405
406 double oneOverDet = 1.0/det;
407 double mn1OverDet = - oneOverDet;
408
409 m[M00] = Det4_1234_1234 * oneOverDet;
410 m[M01] = Det4_0234_1234 * mn1OverDet;
411 m[M02] = Det4_0134_1234 * oneOverDet;
412 m[M03] = Det4_0124_1234 * mn1OverDet;
413 m[M04] = Det4_0123_1234 * oneOverDet;
414
415 m[M10] = Det4_1234_0234 * mn1OverDet;
416 m[M11] = Det4_0234_0234 * oneOverDet;
417 m[M12] = Det4_0134_0234 * mn1OverDet;
418 m[M13] = Det4_0124_0234 * oneOverDet;
419 m[M14] = Det4_0123_0234 * mn1OverDet;
420
421 m[M20] = Det4_1234_0134 * oneOverDet;
422 m[M21] = Det4_0234_0134 * mn1OverDet;
423 m[M22] = Det4_0134_0134 * oneOverDet;
424 m[M23] = Det4_0124_0134 * mn1OverDet;
425 m[M24] = Det4_0123_0134 * oneOverDet;
426
427 m[M30] = Det4_1234_0124 * mn1OverDet;
428 m[M31] = Det4_0234_0124 * oneOverDet;
429 m[M32] = Det4_0134_0124 * mn1OverDet;
430 m[M33] = Det4_0124_0124 * oneOverDet;
431 m[M34] = Det4_0123_0124 * mn1OverDet;
432
433 m[M40] = Det4_1234_0123 * oneOverDet;
434 m[M41] = Det4_0234_0123 * mn1OverDet;
435 m[M42] = Det4_0134_0123 * oneOverDet;
436 m[M43] = Det4_0124_0123 * mn1OverDet;
437 m[M44] = Det4_0123_0123 * oneOverDet;
438
439 return;
440}
441␌
442void HepMatrix::invertHaywood6 (int & ifail) {
443
444 ifail = 0;
445
446 // Find all NECESSARY 2x2 dets: (45 of them)
447
448 double Det2_34_01 = m[A30]*m[A41] - m[A31]*m[A40];
449 double Det2_34_02 = m[A30]*m[A42] - m[A32]*m[A40];
450 double Det2_34_03 = m[A30]*m[A43] - m[A33]*m[A40];
451 double Det2_34_04 = m[A30]*m[A44] - m[A34]*m[A40];
452 double Det2_34_05 = m[A30]*m[A45] - m[A35]*m[A40]; //
453 double Det2_34_12 = m[A31]*m[A42] - m[A32]*m[A41];
454 double Det2_34_13 = m[A31]*m[A43] - m[A33]*m[A41];
455 double Det2_34_14 = m[A31]*m[A44] - m[A34]*m[A41];
456 double Det2_34_15 = m[A31]*m[A45] - m[A35]*m[A41]; //
457 double Det2_34_23 = m[A32]*m[A43] - m[A33]*m[A42];
458 double Det2_34_24 = m[A32]*m[A44] - m[A34]*m[A42];
459 double Det2_34_25 = m[A32]*m[A45] - m[A35]*m[A42]; //
460 double Det2_34_34 = m[A33]*m[A44] - m[A34]*m[A43];
461 double Det2_34_35 = m[A33]*m[A45] - m[A35]*m[A43]; //
462 double Det2_34_45 = m[A34]*m[A45] - m[A35]*m[A44]; //
463 double Det2_35_01 = m[A30]*m[A51] - m[A31]*m[A50];
464 double Det2_35_02 = m[A30]*m[A52] - m[A32]*m[A50];
465 double Det2_35_03 = m[A30]*m[A53] - m[A33]*m[A50];
466 double Det2_35_04 = m[A30]*m[A54] - m[A34]*m[A50];
467 double Det2_35_05 = m[A30]*m[A55] - m[A35]*m[A50];
468 double Det2_35_12 = m[A31]*m[A52] - m[A32]*m[A51];
469 double Det2_35_13 = m[A31]*m[A53] - m[A33]*m[A51];
470 double Det2_35_14 = m[A31]*m[A54] - m[A34]*m[A51];
471 double Det2_35_15 = m[A31]*m[A55] - m[A35]*m[A51];
472 double Det2_35_23 = m[A32]*m[A53] - m[A33]*m[A52];
473 double Det2_35_24 = m[A32]*m[A54] - m[A34]*m[A52];
474 double Det2_35_25 = m[A32]*m[A55] - m[A35]*m[A52];
475 double Det2_35_34 = m[A33]*m[A54] - m[A34]*m[A53];
476 double Det2_35_35 = m[A33]*m[A55] - m[A35]*m[A53];
477 double Det2_35_45 = m[A34]*m[A55] - m[A35]*m[A54]; //
478 double Det2_45_01 = m[A40]*m[A51] - m[A41]*m[A50];
479 double Det2_45_02 = m[A40]*m[A52] - m[A42]*m[A50];
480 double Det2_45_03 = m[A40]*m[A53] - m[A43]*m[A50];
481 double Det2_45_04 = m[A40]*m[A54] - m[A44]*m[A50];
482 double Det2_45_05 = m[A40]*m[A55] - m[A45]*m[A50];
483 double Det2_45_12 = m[A41]*m[A52] - m[A42]*m[A51];
484 double Det2_45_13 = m[A41]*m[A53] - m[A43]*m[A51];
485 double Det2_45_14 = m[A41]*m[A54] - m[A44]*m[A51];
486 double Det2_45_15 = m[A41]*m[A55] - m[A45]*m[A51];
487 double Det2_45_23 = m[A42]*m[A53] - m[A43]*m[A52];
488 double Det2_45_24 = m[A42]*m[A54] - m[A44]*m[A52];
489 double Det2_45_25 = m[A42]*m[A55] - m[A45]*m[A52];
490 double Det2_45_34 = m[A43]*m[A54] - m[A44]*m[A53];
491 double Det2_45_35 = m[A43]*m[A55] - m[A45]*m[A53];
492 double Det2_45_45 = m[A44]*m[A55] - m[A45]*m[A54];
493
494 // Find all NECESSARY 3x3 dets: (80 of them)
495
496 double Det3_234_012 = m[A20]*Det2_34_12 - m[A21]*Det2_34_02
497 + m[A22]*Det2_34_01;
498 double Det3_234_013 = m[A20]*Det2_34_13 - m[A21]*Det2_34_03
499 + m[A23]*Det2_34_01;
500 double Det3_234_014 = m[A20]*Det2_34_14 - m[A21]*Det2_34_04
501 + m[A24]*Det2_34_01;
502 double Det3_234_015 = m[A20]*Det2_34_15 - m[A21]*Det2_34_05
503 + m[A25]*Det2_34_01; //
504 double Det3_234_023 = m[A20]*Det2_34_23 - m[A22]*Det2_34_03
505 + m[A23]*Det2_34_02;
506 double Det3_234_024 = m[A20]*Det2_34_24 - m[A22]*Det2_34_04
507 + m[A24]*Det2_34_02;
508 double Det3_234_025 = m[A20]*Det2_34_25 - m[A22]*Det2_34_05
509 + m[A25]*Det2_34_02; //
510 double Det3_234_034 = m[A20]*Det2_34_34 - m[A23]*Det2_34_04
511 + m[A24]*Det2_34_03;
512 double Det3_234_035 = m[A20]*Det2_34_35 - m[A23]*Det2_34_05
513 + m[A25]*Det2_34_03; //
514 double Det3_234_045 = m[A20]*Det2_34_45 - m[A24]*Det2_34_05
515 + m[A25]*Det2_34_04; //
516 double Det3_234_123 = m[A21]*Det2_34_23 - m[A22]*Det2_34_13
517 + m[A23]*Det2_34_12;
518 double Det3_234_124 = m[A21]*Det2_34_24 - m[A22]*Det2_34_14
519 + m[A24]*Det2_34_12;
520 double Det3_234_125 = m[A21]*Det2_34_25 - m[A22]*Det2_34_15
521 + m[A25]*Det2_34_12; //
522 double Det3_234_134 = m[A21]*Det2_34_34 - m[A23]*Det2_34_14
523 + m[A24]*Det2_34_13;
524 double Det3_234_135 = m[A21]*Det2_34_35 - m[A23]*Det2_34_15
525 + m[A25]*Det2_34_13; //
526 double Det3_234_145 = m[A21]*Det2_34_45 - m[A24]*Det2_34_15
527 + m[A25]*Det2_34_14; //
528 double Det3_234_234 = m[A22]*Det2_34_34 - m[A23]*Det2_34_24
529 + m[A24]*Det2_34_23;
530 double Det3_234_235 = m[A22]*Det2_34_35 - m[A23]*Det2_34_25
531 + m[A25]*Det2_34_23; //
532 double Det3_234_245 = m[A22]*Det2_34_45 - m[A24]*Det2_34_25
533 + m[A25]*Det2_34_24; //
534 double Det3_234_345 = m[A23]*Det2_34_45 - m[A24]*Det2_34_35
535 + m[A25]*Det2_34_34; //
536 double Det3_235_012 = m[A20]*Det2_35_12 - m[A21]*Det2_35_02
537 + m[A22]*Det2_35_01;
538 double Det3_235_013 = m[A20]*Det2_35_13 - m[A21]*Det2_35_03
539 + m[A23]*Det2_35_01;
540 double Det3_235_014 = m[A20]*Det2_35_14 - m[A21]*Det2_35_04
541 + m[A24]*Det2_35_01;
542 double Det3_235_015 = m[A20]*Det2_35_15 - m[A21]*Det2_35_05
543 + m[A25]*Det2_35_01;
544 double Det3_235_023 = m[A20]*Det2_35_23 - m[A22]*Det2_35_03
545 + m[A23]*Det2_35_02;
546 double Det3_235_024 = m[A20]*Det2_35_24 - m[A22]*Det2_35_04
547 + m[A24]*Det2_35_02;
548 double Det3_235_025 = m[A20]*Det2_35_25 - m[A22]*Det2_35_05
549 + m[A25]*Det2_35_02;
550 double Det3_235_034 = m[A20]*Det2_35_34 - m[A23]*Det2_35_04
551 + m[A24]*Det2_35_03;
552 double Det3_235_035 = m[A20]*Det2_35_35 - m[A23]*Det2_35_05
553 + m[A25]*Det2_35_03;
554 double Det3_235_045 = m[A20]*Det2_35_45 - m[A24]*Det2_35_05
555 + m[A25]*Det2_35_04; //
556 double Det3_235_123 = m[A21]*Det2_35_23 - m[A22]*Det2_35_13
557 + m[A23]*Det2_35_12;
558 double Det3_235_124 = m[A21]*Det2_35_24 - m[A22]*Det2_35_14
559 + m[A24]*Det2_35_12;
560 double Det3_235_125 = m[A21]*Det2_35_25 - m[A22]*Det2_35_15
561 + m[A25]*Det2_35_12;
562 double Det3_235_134 = m[A21]*Det2_35_34 - m[A23]*Det2_35_14
563 + m[A24]*Det2_35_13;
564 double Det3_235_135 = m[A21]*Det2_35_35 - m[A23]*Det2_35_15
565 + m[A25]*Det2_35_13;
566 double Det3_235_145 = m[A21]*Det2_35_45 - m[A24]*Det2_35_15
567 + m[A25]*Det2_35_14; //
568 double Det3_235_234 = m[A22]*Det2_35_34 - m[A23]*Det2_35_24
569 + m[A24]*Det2_35_23;
570 double Det3_235_235 = m[A22]*Det2_35_35 - m[A23]*Det2_35_25
571 + m[A25]*Det2_35_23;
572 double Det3_235_245 = m[A22]*Det2_35_45 - m[A24]*Det2_35_25
573 + m[A25]*Det2_35_24; //
574 double Det3_235_345 = m[A23]*Det2_35_45 - m[A24]*Det2_35_35
575 + m[A25]*Det2_35_34; //
576 double Det3_245_012 = m[A20]*Det2_45_12 - m[A21]*Det2_45_02
577 + m[A22]*Det2_45_01;
578 double Det3_245_013 = m[A20]*Det2_45_13 - m[A21]*Det2_45_03
579 + m[A23]*Det2_45_01;
580 double Det3_245_014 = m[A20]*Det2_45_14 - m[A21]*Det2_45_04
581 + m[A24]*Det2_45_01;
582 double Det3_245_015 = m[A20]*Det2_45_15 - m[A21]*Det2_45_05
583 + m[A25]*Det2_45_01;
584 double Det3_245_023 = m[A20]*Det2_45_23 - m[A22]*Det2_45_03
585 + m[A23]*Det2_45_02;
586 double Det3_245_024 = m[A20]*Det2_45_24 - m[A22]*Det2_45_04
587 + m[A24]*Det2_45_02;
588 double Det3_245_025 = m[A20]*Det2_45_25 - m[A22]*Det2_45_05
589 + m[A25]*Det2_45_02;
590 double Det3_245_034 = m[A20]*Det2_45_34 - m[A23]*Det2_45_04
591 + m[A24]*Det2_45_03;
592 double Det3_245_035 = m[A20]*Det2_45_35 - m[A23]*Det2_45_05
593 + m[A25]*Det2_45_03;
594 double Det3_245_045 = m[A20]*Det2_45_45 - m[A24]*Det2_45_05
595 + m[A25]*Det2_45_04;
596 double Det3_245_123 = m[A21]*Det2_45_23 - m[A22]*Det2_45_13
597 + m[A23]*Det2_45_12;
598 double Det3_245_124 = m[A21]*Det2_45_24 - m[A22]*Det2_45_14
599 + m[A24]*Det2_45_12;
600 double Det3_245_125 = m[A21]*Det2_45_25 - m[A22]*Det2_45_15
601 + m[A25]*Det2_45_12;
602 double Det3_245_134 = m[A21]*Det2_45_34 - m[A23]*Det2_45_14
603 + m[A24]*Det2_45_13;
604 double Det3_245_135 = m[A21]*Det2_45_35 - m[A23]*Det2_45_15
605 + m[A25]*Det2_45_13;
606 double Det3_245_145 = m[A21]*Det2_45_45 - m[A24]*Det2_45_15
607 + m[A25]*Det2_45_14;
608 double Det3_245_234 = m[A22]*Det2_45_34 - m[A23]*Det2_45_24
609 + m[A24]*Det2_45_23;
610 double Det3_245_235 = m[A22]*Det2_45_35 - m[A23]*Det2_45_25
611 + m[A25]*Det2_45_23;
612 double Det3_245_245 = m[A22]*Det2_45_45 - m[A24]*Det2_45_25
613 + m[A25]*Det2_45_24;
614 double Det3_245_345 = m[A23]*Det2_45_45 - m[A24]*Det2_45_35
615 + m[A25]*Det2_45_34; //
616 double Det3_345_012 = m[A30]*Det2_45_12 - m[A31]*Det2_45_02
617 + m[A32]*Det2_45_01;
618 double Det3_345_013 = m[A30]*Det2_45_13 - m[A31]*Det2_45_03
619 + m[A33]*Det2_45_01;
620 double Det3_345_014 = m[A30]*Det2_45_14 - m[A31]*Det2_45_04
621 + m[A34]*Det2_45_01;
622 double Det3_345_015 = m[A30]*Det2_45_15 - m[A31]*Det2_45_05
623 + m[A35]*Det2_45_01;
624 double Det3_345_023 = m[A30]*Det2_45_23 - m[A32]*Det2_45_03
625 + m[A33]*Det2_45_02;
626 double Det3_345_024 = m[A30]*Det2_45_24 - m[A32]*Det2_45_04
627 + m[A34]*Det2_45_02;
628 double Det3_345_025 = m[A30]*Det2_45_25 - m[A32]*Det2_45_05
629 + m[A35]*Det2_45_02;
630 double Det3_345_034 = m[A30]*Det2_45_34 - m[A33]*Det2_45_04
631 + m[A34]*Det2_45_03;
632 double Det3_345_035 = m[A30]*Det2_45_35 - m[A33]*Det2_45_05
633 + m[A35]*Det2_45_03;
634 double Det3_345_045 = m[A30]*Det2_45_45 - m[A34]*Det2_45_05
635 + m[A35]*Det2_45_04;
636 double Det3_345_123 = m[A31]*Det2_45_23 - m[A32]*Det2_45_13
637 + m[A33]*Det2_45_12;
638 double Det3_345_124 = m[A31]*Det2_45_24 - m[A32]*Det2_45_14
639 + m[A34]*Det2_45_12;
640 double Det3_345_125 = m[A31]*Det2_45_25 - m[A32]*Det2_45_15
641 + m[A35]*Det2_45_12;
642 double Det3_345_134 = m[A31]*Det2_45_34 - m[A33]*Det2_45_14
643 + m[A34]*Det2_45_13;
644 double Det3_345_135 = m[A31]*Det2_45_35 - m[A33]*Det2_45_15
645 + m[A35]*Det2_45_13;
646 double Det3_345_145 = m[A31]*Det2_45_45 - m[A34]*Det2_45_15
647 + m[A35]*Det2_45_14;
648 double Det3_345_234 = m[A32]*Det2_45_34 - m[A33]*Det2_45_24
649 + m[A34]*Det2_45_23;
650 double Det3_345_235 = m[A32]*Det2_45_35 - m[A33]*Det2_45_25
651 + m[A35]*Det2_45_23;
652 double Det3_345_245 = m[A32]*Det2_45_45 - m[A34]*Det2_45_25
653 + m[A35]*Det2_45_24;
654 double Det3_345_345 = m[A33]*Det2_45_45 - m[A34]*Det2_45_35
655 + m[A35]*Det2_45_34;
656
657 // Find all NECESSARY 4x4 dets: (75 of them)
658
659 double Det4_1234_0123 = m[A10]*Det3_234_123 - m[A11]*Det3_234_023
660 + m[A12]*Det3_234_013 - m[A13]*Det3_234_012;
661 double Det4_1234_0124 = m[A10]*Det3_234_124 - m[A11]*Det3_234_024
662 + m[A12]*Det3_234_014 - m[A14]*Det3_234_012;
663 double Det4_1234_0125 = m[A10]*Det3_234_125 - m[A11]*Det3_234_025
664 + m[A12]*Det3_234_015 - m[A15]*Det3_234_012; //
665 double Det4_1234_0134 = m[A10]*Det3_234_134 - m[A11]*Det3_234_034
666 + m[A13]*Det3_234_014 - m[A14]*Det3_234_013;
667 double Det4_1234_0135 = m[A10]*Det3_234_135 - m[A11]*Det3_234_035
668 + m[A13]*Det3_234_015 - m[A15]*Det3_234_013; //
669 double Det4_1234_0145 = m[A10]*Det3_234_145 - m[A11]*Det3_234_045
670 + m[A14]*Det3_234_015 - m[A15]*Det3_234_014; //
671 double Det4_1234_0234 = m[A10]*Det3_234_234 - m[A12]*Det3_234_034
672 + m[A13]*Det3_234_024 - m[A14]*Det3_234_023;
673 double Det4_1234_0235 = m[A10]*Det3_234_235 - m[A12]*Det3_234_035
674 + m[A13]*Det3_234_025 - m[A15]*Det3_234_023; //
675 double Det4_1234_0245 = m[A10]*Det3_234_245 - m[A12]*Det3_234_045
676 + m[A14]*Det3_234_025 - m[A15]*Det3_234_024; //
677 double Det4_1234_0345 = m[A10]*Det3_234_345 - m[A13]*Det3_234_045
678 + m[A14]*Det3_234_035 - m[A15]*Det3_234_034; //
679 double Det4_1234_1234 = m[A11]*Det3_234_234 - m[A12]*Det3_234_134
680 + m[A13]*Det3_234_124 - m[A14]*Det3_234_123;
681 double Det4_1234_1235 = m[A11]*Det3_234_235 - m[A12]*Det3_234_135
682 + m[A13]*Det3_234_125 - m[A15]*Det3_234_123; //
683 double Det4_1234_1245 = m[A11]*Det3_234_245 - m[A12]*Det3_234_145
684 + m[A14]*Det3_234_125 - m[A15]*Det3_234_124; //
685 double Det4_1234_1345 = m[A11]*Det3_234_345 - m[A13]*Det3_234_145
686 + m[A14]*Det3_234_135 - m[A15]*Det3_234_134; //
687 double Det4_1234_2345 = m[A12]*Det3_234_345 - m[A13]*Det3_234_245
688 + m[A14]*Det3_234_235 - m[A15]*Det3_234_234; //
689 double Det4_1235_0123 = m[A10]*Det3_235_123 - m[A11]*Det3_235_023
690 + m[A12]*Det3_235_013 - m[A13]*Det3_235_012;
691 double Det4_1235_0124 = m[A10]*Det3_235_124 - m[A11]*Det3_235_024
692 + m[A12]*Det3_235_014 - m[A14]*Det3_235_012;
693 double Det4_1235_0125 = m[A10]*Det3_235_125 - m[A11]*Det3_235_025
694 + m[A12]*Det3_235_015 - m[A15]*Det3_235_012;
695 double Det4_1235_0134 = m[A10]*Det3_235_134 - m[A11]*Det3_235_034
696 + m[A13]*Det3_235_014 - m[A14]*Det3_235_013;
697 double Det4_1235_0135 = m[A10]*Det3_235_135 - m[A11]*Det3_235_035
698 + m[A13]*Det3_235_015 - m[A15]*Det3_235_013;
699 double Det4_1235_0145 = m[A10]*Det3_235_145 - m[A11]*Det3_235_045
700 + m[A14]*Det3_235_015 - m[A15]*Det3_235_014; //
701 double Det4_1235_0234 = m[A10]*Det3_235_234 - m[A12]*Det3_235_034
702 + m[A13]*Det3_235_024 - m[A14]*Det3_235_023;
703 double Det4_1235_0235 = m[A10]*Det3_235_235 - m[A12]*Det3_235_035
704 + m[A13]*Det3_235_025 - m[A15]*Det3_235_023;
705 double Det4_1235_0245 = m[A10]*Det3_235_245 - m[A12]*Det3_235_045
706 + m[A14]*Det3_235_025 - m[A15]*Det3_235_024; //
707 double Det4_1235_0345 = m[A10]*Det3_235_345 - m[A13]*Det3_235_045
708 + m[A14]*Det3_235_035 - m[A15]*Det3_235_034; //
709 double Det4_1235_1234 = m[A11]*Det3_235_234 - m[A12]*Det3_235_134
710 + m[A13]*Det3_235_124 - m[A14]*Det3_235_123;
711 double Det4_1235_1235 = m[A11]*Det3_235_235 - m[A12]*Det3_235_135
712 + m[A13]*Det3_235_125 - m[A15]*Det3_235_123;
713 double Det4_1235_1245 = m[A11]*Det3_235_245 - m[A12]*Det3_235_145
714 + m[A14]*Det3_235_125 - m[A15]*Det3_235_124; //
715 double Det4_1235_1345 = m[A11]*Det3_235_345 - m[A13]*Det3_235_145
716 + m[A14]*Det3_235_135 - m[A15]*Det3_235_134; //
717 double Det4_1235_2345 = m[A12]*Det3_235_345 - m[A13]*Det3_235_245
718 + m[A14]*Det3_235_235 - m[A15]*Det3_235_234; //
719 double Det4_1245_0123 = m[A10]*Det3_245_123 - m[A11]*Det3_245_023
720 + m[A12]*Det3_245_013 - m[A13]*Det3_245_012;
721 double Det4_1245_0124 = m[A10]*Det3_245_124 - m[A11]*Det3_245_024
722 + m[A12]*Det3_245_014 - m[A14]*Det3_245_012;
723 double Det4_1245_0125 = m[A10]*Det3_245_125 - m[A11]*Det3_245_025
724 + m[A12]*Det3_245_015 - m[A15]*Det3_245_012;
725 double Det4_1245_0134 = m[A10]*Det3_245_134 - m[A11]*Det3_245_034
726 + m[A13]*Det3_245_014 - m[A14]*Det3_245_013;
727 double Det4_1245_0135 = m[A10]*Det3_245_135 - m[A11]*Det3_245_035
728 + m[A13]*Det3_245_015 - m[A15]*Det3_245_013;
729 double Det4_1245_0145 = m[A10]*Det3_245_145 - m[A11]*Det3_245_045
730 + m[A14]*Det3_245_015 - m[A15]*Det3_245_014;
731 double Det4_1245_0234 = m[A10]*Det3_245_234 - m[A12]*Det3_245_034
732 + m[A13]*Det3_245_024 - m[A14]*Det3_245_023;
733 double Det4_1245_0235 = m[A10]*Det3_245_235 - m[A12]*Det3_245_035
734 + m[A13]*Det3_245_025 - m[A15]*Det3_245_023;
735 double Det4_1245_0245 = m[A10]*Det3_245_245 - m[A12]*Det3_245_045
736 + m[A14]*Det3_245_025 - m[A15]*Det3_245_024;
737 double Det4_1245_0345 = m[A10]*Det3_245_345 - m[A13]*Det3_245_045
738 + m[A14]*Det3_245_035 - m[A15]*Det3_245_034; //
739 double Det4_1245_1234 = m[A11]*Det3_245_234 - m[A12]*Det3_245_134
740 + m[A13]*Det3_245_124 - m[A14]*Det3_245_123;
741 double Det4_1245_1235 = m[A11]*Det3_245_235 - m[A12]*Det3_245_135
742 + m[A13]*Det3_245_125 - m[A15]*Det3_245_123;
743 double Det4_1245_1245 = m[A11]*Det3_245_245 - m[A12]*Det3_245_145
744 + m[A14]*Det3_245_125 - m[A15]*Det3_245_124;
745 double Det4_1245_1345 = m[A11]*Det3_245_345 - m[A13]*Det3_245_145
746 + m[A14]*Det3_245_135 - m[A15]*Det3_245_134; //
747 double Det4_1245_2345 = m[A12]*Det3_245_345 - m[A13]*Det3_245_245
748 + m[A14]*Det3_245_235 - m[A15]*Det3_245_234; //
749 double Det4_1345_0123 = m[A10]*Det3_345_123 - m[A11]*Det3_345_023
750 + m[A12]*Det3_345_013 - m[A13]*Det3_345_012;
751 double Det4_1345_0124 = m[A10]*Det3_345_124 - m[A11]*Det3_345_024
752 + m[A12]*Det3_345_014 - m[A14]*Det3_345_012;
753 double Det4_1345_0125 = m[A10]*Det3_345_125 - m[A11]*Det3_345_025
754 + m[A12]*Det3_345_015 - m[A15]*Det3_345_012;
755 double Det4_1345_0134 = m[A10]*Det3_345_134 - m[A11]*Det3_345_034
756 + m[A13]*Det3_345_014 - m[A14]*Det3_345_013;
757 double Det4_1345_0135 = m[A10]*Det3_345_135 - m[A11]*Det3_345_035
758 + m[A13]*Det3_345_015 - m[A15]*Det3_345_013;
759 double Det4_1345_0145 = m[A10]*Det3_345_145 - m[A11]*Det3_345_045
760 + m[A14]*Det3_345_015 - m[A15]*Det3_345_014;
761 double Det4_1345_0234 = m[A10]*Det3_345_234 - m[A12]*Det3_345_034
762 + m[A13]*Det3_345_024 - m[A14]*Det3_345_023;
763 double Det4_1345_0235 = m[A10]*Det3_345_235 - m[A12]*Det3_345_035
764 + m[A13]*Det3_345_025 - m[A15]*Det3_345_023;
765 double Det4_1345_0245 = m[A10]*Det3_345_245 - m[A12]*Det3_345_045
766 + m[A14]*Det3_345_025 - m[A15]*Det3_345_024;
767 double Det4_1345_0345 = m[A10]*Det3_345_345 - m[A13]*Det3_345_045
768 + m[A14]*Det3_345_035 - m[A15]*Det3_345_034;
769 double Det4_1345_1234 = m[A11]*Det3_345_234 - m[A12]*Det3_345_134
770 + m[A13]*Det3_345_124 - m[A14]*Det3_345_123;
771 double Det4_1345_1235 = m[A11]*Det3_345_235 - m[A12]*Det3_345_135
772 + m[A13]*Det3_345_125 - m[A15]*Det3_345_123;
773 double Det4_1345_1245 = m[A11]*Det3_345_245 - m[A12]*Det3_345_145
774 + m[A14]*Det3_345_125 - m[A15]*Det3_345_124;
775 double Det4_1345_1345 = m[A11]*Det3_345_345 - m[A13]*Det3_345_145
776 + m[A14]*Det3_345_135 - m[A15]*Det3_345_134;
777 double Det4_1345_2345 = m[A12]*Det3_345_345 - m[A13]*Det3_345_245
778 + m[A14]*Det3_345_235 - m[A15]*Det3_345_234; //
779 double Det4_2345_0123 = m[A20]*Det3_345_123 - m[A21]*Det3_345_023
780 + m[A22]*Det3_345_013 - m[A23]*Det3_345_012;
781 double Det4_2345_0124 = m[A20]*Det3_345_124 - m[A21]*Det3_345_024
782 + m[A22]*Det3_345_014 - m[A24]*Det3_345_012;
783 double Det4_2345_0125 = m[A20]*Det3_345_125 - m[A21]*Det3_345_025
784 + m[A22]*Det3_345_015 - m[A25]*Det3_345_012;
785 double Det4_2345_0134 = m[A20]*Det3_345_134 - m[A21]*Det3_345_034
786 + m[A23]*Det3_345_014 - m[A24]*Det3_345_013;
787 double Det4_2345_0135 = m[A20]*Det3_345_135 - m[A21]*Det3_345_035
788 + m[A23]*Det3_345_015 - m[A25]*Det3_345_013;
789 double Det4_2345_0145 = m[A20]*Det3_345_145 - m[A21]*Det3_345_045
790 + m[A24]*Det3_345_015 - m[A25]*Det3_345_014;
791 double Det4_2345_0234 = m[A20]*Det3_345_234 - m[A22]*Det3_345_034
792 + m[A23]*Det3_345_024 - m[A24]*Det3_345_023;
793 double Det4_2345_0235 = m[A20]*Det3_345_235 - m[A22]*Det3_345_035
794 + m[A23]*Det3_345_025 - m[A25]*Det3_345_023;
795 double Det4_2345_0245 = m[A20]*Det3_345_245 - m[A22]*Det3_345_045
796 + m[A24]*Det3_345_025 - m[A25]*Det3_345_024;
797 double Det4_2345_0345 = m[A20]*Det3_345_345 - m[A23]*Det3_345_045
798 + m[A24]*Det3_345_035 - m[A25]*Det3_345_034;
799 double Det4_2345_1234 = m[A21]*Det3_345_234 - m[A22]*Det3_345_134
800 + m[A23]*Det3_345_124 - m[A24]*Det3_345_123;
801 double Det4_2345_1235 = m[A21]*Det3_345_235 - m[A22]*Det3_345_135
802 + m[A23]*Det3_345_125 - m[A25]*Det3_345_123;
803 double Det4_2345_1245 = m[A21]*Det3_345_245 - m[A22]*Det3_345_145
804 + m[A24]*Det3_345_125 - m[A25]*Det3_345_124;
805 double Det4_2345_1345 = m[A21]*Det3_345_345 - m[A23]*Det3_345_145
806 + m[A24]*Det3_345_135 - m[A25]*Det3_345_134;
807 double Det4_2345_2345 = m[A22]*Det3_345_345 - m[A23]*Det3_345_245
808 + m[A24]*Det3_345_235 - m[A25]*Det3_345_234;
809
810 // Find all NECESSARY 5x5 dets: (36 of them)
811
812 double Det5_01234_01234 = m[A00]*Det4_1234_1234 - m[A01]*Det4_1234_0234
813 + m[A02]*Det4_1234_0134 - m[A03]*Det4_1234_0124 + m[A04]*Det4_1234_0123;
814 double Det5_01234_01235 = m[A00]*Det4_1234_1235 - m[A01]*Det4_1234_0235
815 + m[A02]*Det4_1234_0135 - m[A03]*Det4_1234_0125 + m[A05]*Det4_1234_0123; //
816 double Det5_01234_01245 = m[A00]*Det4_1234_1245 - m[A01]*Det4_1234_0245
817 + m[A02]*Det4_1234_0145 - m[A04]*Det4_1234_0125 + m[A05]*Det4_1234_0124; //
818 double Det5_01234_01345 = m[A00]*Det4_1234_1345 - m[A01]*Det4_1234_0345
819 + m[A03]*Det4_1234_0145 - m[A04]*Det4_1234_0135 + m[A05]*Det4_1234_0134; //
820 double Det5_01234_02345 = m[A00]*Det4_1234_2345 - m[A02]*Det4_1234_0345
821 + m[A03]*Det4_1234_0245 - m[A04]*Det4_1234_0235 + m[A05]*Det4_1234_0234; //
822 double Det5_01234_12345 = m[A01]*Det4_1234_2345 - m[A02]*Det4_1234_1345
823 + m[A03]*Det4_1234_1245 - m[A04]*Det4_1234_1235 + m[A05]*Det4_1234_1234; //
824 double Det5_01235_01234 = m[A00]*Det4_1235_1234 - m[A01]*Det4_1235_0234
825 + m[A02]*Det4_1235_0134 - m[A03]*Det4_1235_0124 + m[A04]*Det4_1235_0123;
826 double Det5_01235_01235 = m[A00]*Det4_1235_1235 - m[A01]*Det4_1235_0235
827 + m[A02]*Det4_1235_0135 - m[A03]*Det4_1235_0125 + m[A05]*Det4_1235_0123;
828 double Det5_01235_01245 = m[A00]*Det4_1235_1245 - m[A01]*Det4_1235_0245
829 + m[A02]*Det4_1235_0145 - m[A04]*Det4_1235_0125 + m[A05]*Det4_1235_0124; //
830 double Det5_01235_01345 = m[A00]*Det4_1235_1345 - m[A01]*Det4_1235_0345
831 + m[A03]*Det4_1235_0145 - m[A04]*Det4_1235_0135 + m[A05]*Det4_1235_0134; //
832 double Det5_01235_02345 = m[A00]*Det4_1235_2345 - m[A02]*Det4_1235_0345
833 + m[A03]*Det4_1235_0245 - m[A04]*Det4_1235_0235 + m[A05]*Det4_1235_0234; //
834 double Det5_01235_12345 = m[A01]*Det4_1235_2345 - m[A02]*Det4_1235_1345
835 + m[A03]*Det4_1235_1245 - m[A04]*Det4_1235_1235 + m[A05]*Det4_1235_1234; //
836 double Det5_01245_01234 = m[A00]*Det4_1245_1234 - m[A01]*Det4_1245_0234
837 + m[A02]*Det4_1245_0134 - m[A03]*Det4_1245_0124 + m[A04]*Det4_1245_0123;
838 double Det5_01245_01235 = m[A00]*Det4_1245_1235 - m[A01]*Det4_1245_0235
839 + m[A02]*Det4_1245_0135 - m[A03]*Det4_1245_0125 + m[A05]*Det4_1245_0123;
840 double Det5_01245_01245 = m[A00]*Det4_1245_1245 - m[A01]*Det4_1245_0245
841 + m[A02]*Det4_1245_0145 - m[A04]*Det4_1245_0125 + m[A05]*Det4_1245_0124;
842 double Det5_01245_01345 = m[A00]*Det4_1245_1345 - m[A01]*Det4_1245_0345
843 + m[A03]*Det4_1245_0145 - m[A04]*Det4_1245_0135 + m[A05]*Det4_1245_0134; //
844 double Det5_01245_02345 = m[A00]*Det4_1245_2345 - m[A02]*Det4_1245_0345
845 + m[A03]*Det4_1245_0245 - m[A04]*Det4_1245_0235 + m[A05]*Det4_1245_0234; //
846 double Det5_01245_12345 = m[A01]*Det4_1245_2345 - m[A02]*Det4_1245_1345
847 + m[A03]*Det4_1245_1245 - m[A04]*Det4_1245_1235 + m[A05]*Det4_1245_1234; //
848 double Det5_01345_01234 = m[A00]*Det4_1345_1234 - m[A01]*Det4_1345_0234
849 + m[A02]*Det4_1345_0134 - m[A03]*Det4_1345_0124 + m[A04]*Det4_1345_0123;
850 double Det5_01345_01235 = m[A00]*Det4_1345_1235 - m[A01]*Det4_1345_0235
851 + m[A02]*Det4_1345_0135 - m[A03]*Det4_1345_0125 + m[A05]*Det4_1345_0123;
852 double Det5_01345_01245 = m[A00]*Det4_1345_1245 - m[A01]*Det4_1345_0245
853 + m[A02]*Det4_1345_0145 - m[A04]*Det4_1345_0125 + m[A05]*Det4_1345_0124;
854 double Det5_01345_01345 = m[A00]*Det4_1345_1345 - m[A01]*Det4_1345_0345
855 + m[A03]*Det4_1345_0145 - m[A04]*Det4_1345_0135 + m[A05]*Det4_1345_0134;
856 double Det5_01345_02345 = m[A00]*Det4_1345_2345 - m[A02]*Det4_1345_0345
857 + m[A03]*Det4_1345_0245 - m[A04]*Det4_1345_0235 + m[A05]*Det4_1345_0234; //
858 double Det5_01345_12345 = m[A01]*Det4_1345_2345 - m[A02]*Det4_1345_1345
859 + m[A03]*Det4_1345_1245 - m[A04]*Det4_1345_1235 + m[A05]*Det4_1345_1234; //
860 double Det5_02345_01234 = m[A00]*Det4_2345_1234 - m[A01]*Det4_2345_0234
861 + m[A02]*Det4_2345_0134 - m[A03]*Det4_2345_0124 + m[A04]*Det4_2345_0123;
862 double Det5_02345_01235 = m[A00]*Det4_2345_1235 - m[A01]*Det4_2345_0235
863 + m[A02]*Det4_2345_0135 - m[A03]*Det4_2345_0125 + m[A05]*Det4_2345_0123;
864 double Det5_02345_01245 = m[A00]*Det4_2345_1245 - m[A01]*Det4_2345_0245
865 + m[A02]*Det4_2345_0145 - m[A04]*Det4_2345_0125 + m[A05]*Det4_2345_0124;
866 double Det5_02345_01345 = m[A00]*Det4_2345_1345 - m[A01]*Det4_2345_0345
867 + m[A03]*Det4_2345_0145 - m[A04]*Det4_2345_0135 + m[A05]*Det4_2345_0134;
868 double Det5_02345_02345 = m[A00]*Det4_2345_2345 - m[A02]*Det4_2345_0345
869 + m[A03]*Det4_2345_0245 - m[A04]*Det4_2345_0235 + m[A05]*Det4_2345_0234;
870 double Det5_02345_12345 = m[A01]*Det4_2345_2345 - m[A02]*Det4_2345_1345
871 + m[A03]*Det4_2345_1245 - m[A04]*Det4_2345_1235 + m[A05]*Det4_2345_1234; //
872 double Det5_12345_01234 = m[A10]*Det4_2345_1234 - m[A11]*Det4_2345_0234
873 + m[A12]*Det4_2345_0134 - m[A13]*Det4_2345_0124 + m[A14]*Det4_2345_0123;
874 double Det5_12345_01235 = m[A10]*Det4_2345_1235 - m[A11]*Det4_2345_0235
875 + m[A12]*Det4_2345_0135 - m[A13]*Det4_2345_0125 + m[A15]*Det4_2345_0123;
876 double Det5_12345_01245 = m[A10]*Det4_2345_1245 - m[A11]*Det4_2345_0245
877 + m[A12]*Det4_2345_0145 - m[A14]*Det4_2345_0125 + m[A15]*Det4_2345_0124;
878 double Det5_12345_01345 = m[A10]*Det4_2345_1345 - m[A11]*Det4_2345_0345
879 + m[A13]*Det4_2345_0145 - m[A14]*Det4_2345_0135 + m[A15]*Det4_2345_0134;
880 double Det5_12345_02345 = m[A10]*Det4_2345_2345 - m[A12]*Det4_2345_0345
881 + m[A13]*Det4_2345_0245 - m[A14]*Det4_2345_0235 + m[A15]*Det4_2345_0234;
882 double Det5_12345_12345 = m[A11]*Det4_2345_2345 - m[A12]*Det4_2345_1345
883 + m[A13]*Det4_2345_1245 - m[A14]*Det4_2345_1235 + m[A15]*Det4_2345_1234;
884
885 // Find the determinant
886
887 double det = m[A00]*Det5_12345_12345
888 - m[A01]*Det5_12345_02345
889 + m[A02]*Det5_12345_01345
890 - m[A03]*Det5_12345_01245
891 + m[A04]*Det5_12345_01235
892 - m[A05]*Det5_12345_01234;
893
894 if ( det == 0 ) {
895#ifdef SINGULAR_DIAGNOSTICS
896 std::cerr << "Kramer's rule inversion of a singular 6x6 matrix: "
897 << *this << "\n";
898#endif
899 ifail = 1;
900 return;
901 }
902
903 double oneOverDet = 1.0/det;
904 double mn1OverDet = - oneOverDet;
905
906 m[A00] = Det5_12345_12345*oneOverDet;
907 m[A01] = Det5_02345_12345*mn1OverDet;
908 m[A02] = Det5_01345_12345*oneOverDet;
909 m[A03] = Det5_01245_12345*mn1OverDet;
910 m[A04] = Det5_01235_12345*oneOverDet;
911 m[A05] = Det5_01234_12345*mn1OverDet;
912
913 m[A10] = Det5_12345_02345*mn1OverDet;
914 m[A11] = Det5_02345_02345*oneOverDet;
915 m[A12] = Det5_01345_02345*mn1OverDet;
916 m[A13] = Det5_01245_02345*oneOverDet;
917 m[A14] = Det5_01235_02345*mn1OverDet;
918 m[A15] = Det5_01234_02345*oneOverDet;
919
920 m[A20] = Det5_12345_01345*oneOverDet;
921 m[A21] = Det5_02345_01345*mn1OverDet;
922 m[A22] = Det5_01345_01345*oneOverDet;
923 m[A23] = Det5_01245_01345*mn1OverDet;
924 m[A24] = Det5_01235_01345*oneOverDet;
925 m[A25] = Det5_01234_01345*mn1OverDet;
926
927 m[A30] = Det5_12345_01245*mn1OverDet;
928 m[A31] = Det5_02345_01245*oneOverDet;
929 m[A32] = Det5_01345_01245*mn1OverDet;
930 m[A33] = Det5_01245_01245*oneOverDet;
931 m[A34] = Det5_01235_01245*mn1OverDet;
932 m[A35] = Det5_01234_01245*oneOverDet;
933
934 m[A40] = Det5_12345_01235*oneOverDet;
935 m[A41] = Det5_02345_01235*mn1OverDet;
936 m[A42] = Det5_01345_01235*oneOverDet;
937 m[A43] = Det5_01245_01235*mn1OverDet;
938 m[A44] = Det5_01235_01235*oneOverDet;
939 m[A45] = Det5_01234_01235*mn1OverDet;
940
941 m[A50] = Det5_12345_01234*mn1OverDet;
942 m[A51] = Det5_02345_01234*oneOverDet;
943 m[A52] = Det5_01345_01234*mn1OverDet;
944 m[A53] = Det5_01245_01234*oneOverDet;
945 m[A54] = Det5_01235_01234*mn1OverDet;
946 m[A55] = Det5_01234_01234*oneOverDet;
947
948 return;
949}
950
951
952} // namespace CLHEP
A41
#define A41
Definition: MatrixInvert.cc:48
F01
#define F01
Definition: MatrixInvert.cc:92
F22
#define F22
Definition: MatrixInvert.cc:103
F00
#define F00
Definition: MatrixInvert.cc:91
A20
#define A20
Definition: MatrixInvert.cc:33
A01
#define A01
Definition: MatrixInvert.cc:20
A53
#define A53
Definition: MatrixInvert.cc:57
M00
#define M00
Definition: MatrixInvert.cc:61
F02
#define F02
Definition: MatrixInvert.cc:93
A23
#define A23
Definition: MatrixInvert.cc:36
M34
#define M34
Definition: MatrixInvert.cc:83
A45
#define A45
Definition: MatrixInvert.cc:52
A13
#define A13
Definition: MatrixInvert.cc:29
A00
#define A00
Definition: MatrixInvert.cc:19
M03
#define M03
Definition: MatrixInvert.cc:64
A54
#define A54
Definition: MatrixInvert.cc:58
F31
#define F31
Definition: MatrixInvert.cc:107
A40
#define A40
Definition: MatrixInvert.cc:47
F32
#define F32
Definition: MatrixInvert.cc:108
M10
#define M10
Definition: MatrixInvert.cc:67
F21
#define F21
Definition: MatrixInvert.cc:102
F03
#define F03
Definition: MatrixInvert.cc:94
M41
#define M41
Definition: MatrixInvert.cc:86
M20
#define M20
Definition: MatrixInvert.cc:73
M42
#define M42
Definition: MatrixInvert.cc:87
A25
#define A25
Definition: MatrixInvert.cc:38
M13
#define M13
Definition: MatrixInvert.cc:70
A02
#define A02
Definition: MatrixInvert.cc:21
F30
#define F30
Definition: MatrixInvert.cc:106
A24
#define A24
Definition: MatrixInvert.cc:37
M33
#define M33
Definition: MatrixInvert.cc:82
A22
#define A22
Definition: MatrixInvert.cc:35
M04
#define M04
Definition: MatrixInvert.cc:65
A04
#define A04
Definition: MatrixInvert.cc:23
A30
#define A30
Definition: MatrixInvert.cc:40
M23
#define M23
Definition: MatrixInvert.cc:76
A12
#define A12
Definition: MatrixInvert.cc:28
A55
#define A55
Definition: MatrixInvert.cc:59
M11
#define M11
Definition: MatrixInvert.cc:68
A35
#define A35
Definition: MatrixInvert.cc:45
M44
#define M44
Definition: MatrixInvert.cc:89
A50
#define A50
Definition: MatrixInvert.cc:54
A03
#define A03
Definition: MatrixInvert.cc:22
A14
#define A14
Definition: MatrixInvert.cc:30
F11
#define F11
Definition: MatrixInvert.cc:97
M31
#define M31
Definition: MatrixInvert.cc:80
F10
#define F10
Definition: MatrixInvert.cc:96
A51
#define A51
Definition: MatrixInvert.cc:55
M02
#define M02
Definition: MatrixInvert.cc:63
M21
#define M21
Definition: MatrixInvert.cc:74
A21
#define A21
Definition: MatrixInvert.cc:34
M01
#define M01
Definition: MatrixInvert.cc:62
A11
#define A11
Definition: MatrixInvert.cc:27
A10
#define A10
Definition: MatrixInvert.cc:26
A44
#define A44
Definition: MatrixInvert.cc:51
A05
#define A05
Definition: MatrixInvert.cc:24
A32
#define A32
Definition: MatrixInvert.cc:42
A31
#define A31
Definition: MatrixInvert.cc:41
F33
#define F33
Definition: MatrixInvert.cc:109
M24
#define M24
Definition: MatrixInvert.cc:77
F20
#define F20
Definition: MatrixInvert.cc:101
A33
#define A33
Definition: MatrixInvert.cc:43
M12
#define M12
Definition: MatrixInvert.cc:69
M22
#define M22
Definition: MatrixInvert.cc:75
A42
#define A42
Definition: MatrixInvert.cc:49
F12
#define F12
Definition: MatrixInvert.cc:98
A52
#define A52
Definition: MatrixInvert.cc:56
M43
#define M43
Definition: MatrixInvert.cc:88
F23
#define F23
Definition: MatrixInvert.cc:104
A15
#define A15
Definition: MatrixInvert.cc:31
A34
#define A34
Definition: MatrixInvert.cc:44
M14
#define M14
Definition: MatrixInvert.cc:71
M30
#define M30
Definition: MatrixInvert.cc:79
M40
#define M40
Definition: MatrixInvert.cc:85
A43
#define A43
Definition: MatrixInvert.cc:50
M32
#define M32
Definition: MatrixInvert.cc:81
F13
#define F13
Definition: MatrixInvert.cc:99
CLHEP::HepMatrix::invertHaywood4
virtual void invertHaywood4(int &ierr)
Definition: MatrixInvert.cc:112
CLHEP::HepMatrix::invertHaywood6
virtual void invertHaywood6(int &ierr)
Definition: MatrixInvert.cc:442
CLHEP::HepMatrix::invertHaywood5
virtual void invertHaywood5(int &ierr)
Definition: MatrixInvert.cc:216