Geant4-11
Data Structures | Public Member Functions | Static Public Member Functions | Protected Member Functions | Private Member Functions | Private Attributes | Friends
G4ErrorMatrix Class Reference

#include <G4ErrorMatrix.hh>

Data Structures

class  G4ErrorMatrix_row
 
class  G4ErrorMatrix_row_const
 

Public Member Functions

G4ErrorMatrix apply (G4double(*f)(G4double, G4int, G4int)) const
 
G4double determinant () const
 
 G4ErrorMatrix ()
 
 G4ErrorMatrix (const G4ErrorMatrix &m1)
 
 G4ErrorMatrix (const G4ErrorSymMatrix &m1)
 
 G4ErrorMatrix (G4ErrorMatrix &&)=default
 
 G4ErrorMatrix (G4int p, G4int q)
 
 G4ErrorMatrix (G4int p, G4int q, G4int i)
 
G4ErrorMatrix inverse (G4int &ierr) const
 
virtual void invert (G4int &ierr)
 
virtual G4int num_col () const
 
virtual G4int num_row () const
 
virtual G4double & operator() (G4int row, G4int col)
 
virtual const G4double & operator() (G4int row, G4int col) const
 
G4ErrorMatrix & operator*= (G4double t)
 
G4ErrorMatrix & operator+= (const G4ErrorMatrix &m2)
 
G4ErrorMatrix & operator+= (const G4ErrorSymMatrix &m2)
 
G4ErrorMatrix operator- () const
 
G4ErrorMatrix & operator-= (const G4ErrorMatrix &m2)
 
G4ErrorMatrix & operator-= (const G4ErrorSymMatrix &m2)
 
G4ErrorMatrix & operator/= (G4double t)
 
G4ErrorMatrix & operator= (const G4ErrorMatrix &m2)
 
G4ErrorMatrix & operator= (const G4ErrorSymMatrix &m2)
 
G4ErrorMatrix & operator= (G4ErrorMatrix &&)=default
 
G4ErrorMatrix_row operator[] (G4int)
 
const G4ErrorMatrix_row_const operator[] (G4int) const
 
G4ErrorMatrix sub (G4int min_row, G4int max_row, G4int min_col, G4int max_col) const
 
void sub (G4int row, G4int col, const G4ErrorMatrix &m1)
 
G4ErrorMatrix T () const
 
G4double trace () const
 
virtual ~G4ErrorMatrix ()
 

Static Public Member Functions

static void error (const char *s)
 

Protected Member Functions

virtual void invertHaywood4 (G4int &ierr)
 
virtual void invertHaywood5 (G4int &ierr)
 
virtual void invertHaywood6 (G4int &ierr)
 
virtual G4int num_size () const
 

Private Member Functions

G4int dfact_matrix (G4double &det, G4int *ir)
 
G4int dfinv_matrix (G4int *ir)
 

Private Attributes

std::vector< G4double > m
 
G4int ncol
 
G4int nrow
 
G4int size
 

Friends

void back_solve (const G4ErrorMatrix &R, G4ErrorMatrix *b)
 
void col_givens (G4ErrorMatrix *A, G4double c, G4double s, G4int k1, G4int k2, G4int rowmin, G4int rowmax)
 
void col_house (G4ErrorMatrix *, const G4ErrorMatrix &, G4double, G4int, G4int, G4int, G4int)
 
class G4ErrorMatrix_row
 
class G4ErrorMatrix_row_const
 
class G4ErrorSymMatrix
 
void house_with_update (G4ErrorMatrix *a, G4ErrorMatrix *v, G4int row, G4int col)
 
void house_with_update (G4ErrorMatrix *a, G4int row, G4int col)
 
void house_with_update2 (G4ErrorSymMatrix *a, G4ErrorMatrix *v, G4int row, G4int col)
 
G4ErrorMatrix operator* (const G4ErrorMatrix &m1, const G4ErrorMatrix &m2)
 
G4ErrorMatrix operator* (const G4ErrorMatrix &m1, const G4ErrorSymMatrix &m2)
 
G4ErrorMatrix operator* (const G4ErrorSymMatrix &m1, const G4ErrorMatrix &m2)
 
G4ErrorMatrix operator* (const G4ErrorSymMatrix &m1, const G4ErrorSymMatrix &m2)
 
G4ErrorMatrix operator+ (const G4ErrorMatrix &m1, const G4ErrorMatrix &m2)
 
G4ErrorMatrix operator- (const G4ErrorMatrix &m1, const G4ErrorMatrix &m2)
 
G4ErrorMatrix qr_solve (G4ErrorMatrix *, const G4ErrorMatrix &b)
 
void row_givens (G4ErrorMatrix *A, G4double c, G4double s, G4int k1, G4int k2, G4int colmin, G4int colmax)
 
void row_house (G4ErrorMatrix *, const G4ErrorMatrix &, G4double, G4int, G4int, G4int, G4int)
 
void tridiagonal (G4ErrorSymMatrix *a, G4ErrorMatrix *hsm)
 

Detailed Description

Definition at line 46 of file G4ErrorMatrix.hh.

Constructor & Destructor Documentation

◆ G4ErrorMatrix() [1/6]

G4ErrorMatrix::G4ErrorMatrix ( )

◆ G4ErrorMatrix() [2/6]

G4ErrorMatrix::G4ErrorMatrix ( G4int  p,
G4int  q 
)

Definition at line 78 of file G4ErrorMatrix.cc.

79 : m(p * q)
80 , nrow(p)
81 , ncol(q)
82{
83 size = nrow * ncol;
84}
G4ErrorMatrix::m
std::vector< G4double > m
Definition: G4ErrorMatrix.hh:222

References ncol, nrow, and size.

◆ G4ErrorMatrix() [3/6]

G4ErrorMatrix::G4ErrorMatrix ( G4int  p,
G4int  q,
G4int  i 
)

Definition at line 86 of file G4ErrorMatrix.cc.

87 : m(p * q)
88 , nrow(p)
89 , ncol(q)
90{
91 size = nrow * ncol;
92
93 if(size > 0)
94 {
95 switch(init)
96 {
97 case 0:
98 break;
99
100 case 1:
101 {
102 if(ncol == nrow)
103 {
104 G4ErrorMatrixIter a = m.begin();
105 G4ErrorMatrixIter b = m.end();
106 for(; a < b; a += (ncol + 1))
107 *a = 1.0;
108 }
109 else
110 {
111 error("Invalid dimension in G4ErrorMatrix(G4int,G4int,1).");
112 }
113 break;
114 }
115 default:
116 error("G4ErrorMatrix: initialization must be either 0 or 1.");
117 }
118 }
119}
G4ErrorMatrixIter
std::vector< G4double >::iterator G4ErrorMatrixIter
Definition: G4ErrorMatrix.hh:41
G4ErrorMatrix::error
static void error(const char *s)
Definition: G4ErrorMatrix.cc:855

References error(), m, ncol, nrow, and size.

◆ G4ErrorMatrix() [4/6]

G4ErrorMatrix::G4ErrorMatrix ( const G4ErrorMatrix &  m1)

Definition at line 126 of file G4ErrorMatrix.cc.

127 : m(mat1.size)
128 , nrow(mat1.nrow)
129 , ncol(mat1.ncol)
130 , size(mat1.size)
131{
132 m = mat1.m;
133}

References m.

◆ G4ErrorMatrix() [5/6]

G4ErrorMatrix::G4ErrorMatrix ( G4ErrorMatrix &&  )
default

◆ G4ErrorMatrix() [6/6]

G4ErrorMatrix::G4ErrorMatrix ( const G4ErrorSymMatrix &  m1)

Definition at line 135 of file G4ErrorMatrix.cc.

136 : m(mat1.nrow * mat1.nrow)
137 , nrow(mat1.nrow)
138 , ncol(mat1.nrow)
139{
140 size = nrow * ncol;
141
142 G4int n = ncol;
143 G4ErrorMatrixConstIter sjk = mat1.m.begin();
144 G4ErrorMatrixIter m1j = m.begin();
145 G4ErrorMatrixIter mj = m.begin();
146 // j >= k
147 for(G4int j = 1; j <= nrow; j++)
148 {
149 G4ErrorMatrixIter mjk = mj;
150 G4ErrorMatrixIter mkj = m1j;
151 for(G4int k = 1; k <= j; k++)
152 {
153 *(mjk++) = *sjk;
154 if(j != k)
155 *mkj = *sjk;
156 sjk++;
157 mkj += n;
158 }
159 mj += n;
160 m1j++;
161 }
162}
G4ErrorMatrixConstIter
std::vector< G4double >::const_iterator G4ErrorMatrixConstIter
Definition: G4ErrorMatrix.hh:44
G4int
int G4int
Definition: G4Types.hh:85

References m, G4ErrorSymMatrix::m, CLHEP::detail::n, ncol, nrow, and size.

◆ ~G4ErrorMatrix()

G4ErrorMatrix::~G4ErrorMatrix ( )
virtual

Definition at line 124 of file G4ErrorMatrix.cc.

124{}

Member Function Documentation

◆ apply()

G4ErrorMatrix G4ErrorMatrix::apply ( G4double(*)(G4double, G4int, G4int)  f) const

Definition at line 417 of file G4ErrorMatrix.cc.

418{
419 G4ErrorMatrix mret(num_row(), num_col());
420 G4ErrorMatrixConstIter a = m.begin();
421 G4ErrorMatrixIter b = mret.m.begin();
422 for(G4int ir = 1; ir <= num_row(); ir++)
423 {
424 for(G4int ic = 1; ic <= num_col(); ic++)
425 {
426 *(b++) = (*f)(*(a++), ir, ic);
427 }
428 }
429 return mret;
430}
G4ErrorMatrix::num_col
virtual G4int num_col() const
G4ErrorMatrix::num_row
virtual G4int num_row() const

References m, num_col(), and num_row().

◆ determinant()

G4double G4ErrorMatrix::determinant ( ) const

Definition at line 821 of file G4ErrorMatrix.cc.

822{
823 static G4ThreadLocal G4int max_array = 20;
824 static G4ThreadLocal G4int* ir = 0;
825 if(!ir)
826 ir = new G4int[max_array + 1];
827 if(ncol != nrow)
828 {
829 error("G4ErrorMatrix::determinant: G4ErrorMatrix is not NxN");
830 }
831 if(ncol > max_array)
832 {
833 delete[] ir;
834 max_array = nrow;
835 ir = new G4int[max_array + 1];
836 }
837 G4double det;
838 G4ErrorMatrix mt(*this);
839 G4int i = mt.dfact_matrix(det, ir);
840 if(i == 0)
841 return det;
842 return 0;
843}
G4double
double G4double
Definition: G4Types.hh:83
G4ThreadLocal
#define G4ThreadLocal
Definition: tls.hh:77

References dfact_matrix(), error(), G4ThreadLocal, ncol, and nrow.

◆ dfact_matrix()

G4int G4ErrorMatrix::dfact_matrix ( G4double &  det,
G4int *  ir 
)
private

Definition at line 552 of file G4ErrorMatrix.cc.

553{
554 if(ncol != nrow)
555 error("dfact_matrix: G4ErrorMatrix is not NxN");
556
557 G4int ifail, jfail;
558 G4int n = ncol;
559
560 G4double tf;
561 G4double g1 = 1.0e-19, g2 = 1.0e19;
562
563 G4double p, q, t;
564 G4double s11, s12;
565
566 G4double epsilon = 8 * DBL_EPSILON;
567 // could be set to zero (like it was before)
568 // but then the algorithm often doesn't detect
569 // that a matrix is singular
570
571 G4int normal = 0, imposs = -1;
572 G4int jrange = 0, jover = 1, junder = -1;
573 ifail = normal;
574 jfail = jrange;
575 G4int nxch = 0;
576 det = 1.0;
577 G4ErrorMatrixIter mj = m.begin();
578 G4ErrorMatrixIter mjj = mj;
579 for(G4int j = 1; j <= n; j++)
580 {
581 G4int k = j;
582 p = (std::fabs(*mjj));
583 if(j != n)
584 {
585 G4ErrorMatrixIter mij = mj + n + j - 1;
586 for(G4int i = j + 1; i <= n; i++)
587 {
588 q = (std::fabs(*(mij)));
589 if(q > p)
590 {
591 k = i;
592 p = q;
593 }
594 mij += n;
595 }
596 if(k == j)
597 {
598 if(p <= epsilon)
599 {
600 det = 0;
601 ifail = imposs;
602 jfail = jrange;
603 return ifail;
604 }
605 det = -det; // in this case the sign of the determinant
606 // must not change. So I change it twice.
607 }
608 G4ErrorMatrixIter mjl = mj;
609 G4ErrorMatrixIter mkl = m.begin() + (k - 1) * n;
610 for(G4int l = 1; l <= n; l++)
611 {
612 tf = *mjl;
613 *(mjl++) = *mkl;
614 *(mkl++) = tf;
615 }
616 nxch = nxch + 1; // this makes the determinant change its sign
617 ir[nxch] = (((j) << 12) + (k));
618 }
619 else
620 {
621 if(p <= epsilon)
622 {
623 det = 0.0;
624 ifail = imposs;
625 jfail = jrange;
626 return ifail;
627 }
628 }
629 det *= *mjj;
630 *mjj = 1.0 / *mjj;
631 t = (std::fabs(det));
632 if(t < g1)
633 {
634 det = 0.0;
635 if(jfail == jrange)
636 jfail = junder;
637 }
638 else if(t > g2)
639 {
640 det = 1.0;
641 if(jfail == jrange)
642 jfail = jover;
643 }
644 if(j != n)
645 {
646 G4ErrorMatrixIter mk = mj + n;
647 G4ErrorMatrixIter mkjp = mk + j;
648 G4ErrorMatrixIter mjk = mj + j;
649 for(k = j + 1; k <= n; k++)
650 {
651 s11 = -(*mjk);
652 s12 = -(*mkjp);
653 if(j != 1)
654 {
655 G4ErrorMatrixIter mik = m.begin() + k - 1;
656 G4ErrorMatrixIter mijp = m.begin() + j;
657 G4ErrorMatrixIter mki = mk;
658 G4ErrorMatrixIter mji = mj;
659 for(G4int i = 1; i < j; i++)
660 {
661 s11 += (*mik) * (*(mji++));
662 s12 += (*mijp) * (*(mki++));
663 mik += n;
664 mijp += n;
665 }
666 }
667 *(mjk++) = -s11 * (*mjj);
668 *(mkjp) = -(((*(mjj + 1))) * ((*(mkjp - 1))) + (s12));
669 mk += n;
670 mkjp += n;
671 }
672 }
673 mj += n;
674 mjj += (n + 1);
675 }
676 if(nxch % 2 == 1)
677 det = -det;
678 if(jfail != jrange)
679 det = 0.0;
680 ir[n] = nxch;
681 return 0;
682}
epsilon
G4double epsilon(G4double density, G4double temperature)
Definition: G4DNAElectronHoleRecombination.cc:95
CLHEP::normal
static double normal(HepRandomEngine *eptr)
Definition: RandPoisson.cc:79
DBL_EPSILON
#define DBL_EPSILON
Definition: templates.hh:66

References DBL_EPSILON, epsilon(), error(), m, CLHEP::detail::n, ncol, CLHEP::normal(), and nrow.

Referenced by determinant(), G4ErrorSymMatrix::determinant(), and invert().

◆ dfinv_matrix()

G4int G4ErrorMatrix::dfinv_matrix ( G4int *  ir)
private

Definition at line 432 of file G4ErrorMatrix.cc.

433{
434 if(num_col() != num_row())
435 {
436 error("dfinv_matrix: G4ErrorMatrix is not NxN");
437 }
438 G4int n = num_col();
439 if(n == 1)
440 {
441 return 0;
442 }
443
444 G4double s31, s32;
445 G4double s33, s34;
446
447 G4ErrorMatrixIter m11 = m.begin();
448 G4ErrorMatrixIter m12 = m11 + 1;
449 G4ErrorMatrixIter m21 = m11 + n;
450 G4ErrorMatrixIter m22 = m12 + n;
451 *m21 = -(*m22) * (*m11) * (*m21);
452 *m12 = -(*m12);
453 if(n > 2)
454 {
455 G4ErrorMatrixIter mi = m.begin() + 2 * n;
456 G4ErrorMatrixIter mii = m.begin() + 2 * n + 2;
457 G4ErrorMatrixIter mimim = m.begin() + n + 1;
458 for(G4int i = 3; i <= n; i++)
459 {
460 G4int im2 = i - 2;
461 G4ErrorMatrixIter mj = m.begin();
462 G4ErrorMatrixIter mji = mj + i - 1;
463 G4ErrorMatrixIter mij = mi;
464 for(G4int j = 1; j <= im2; j++)
465 {
466 s31 = 0.0;
467 s32 = *mji;
468 G4ErrorMatrixIter mkj = mj + j - 1;
469 G4ErrorMatrixIter mik = mi + j - 1;
470 G4ErrorMatrixIter mjkp = mj + j;
471 G4ErrorMatrixIter mkpi = mj + n + i - 1;
472 for(G4int k = j; k <= im2; k++)
473 {
474 s31 += (*mkj) * (*(mik++));
475 s32 += (*(mjkp++)) * (*mkpi);
476 mkj += n;
477 mkpi += n;
478 }
479 *mij = -(*mii) * (((*(mij - n))) * ((*(mii - 1))) + (s31));
480 *mji = -s32;
481 mj += n;
482 mji += n;
483 mij++;
484 }
485 *(mii - 1) = -(*mii) * (*mimim) * (*(mii - 1));
486 *(mimim + 1) = -(*(mimim + 1));
487 mi += n;
488 mimim += (n + 1);
489 mii += (n + 1);
490 }
491 }
492 G4ErrorMatrixIter mi = m.begin();
493 G4ErrorMatrixIter mii = m.begin();
494 for(G4int i = 1; i < n; i++)
495 {
496 G4int ni = n - i;
497 G4ErrorMatrixIter mij = mi;
498 G4int j;
499 for(j = 1; j <= i; j++)
500 {
501 s33 = *mij;
502 G4ErrorMatrixIter mikj = mi + n + j - 1;
503 G4ErrorMatrixIter miik = mii + 1;
504 G4ErrorMatrixIter min_end = mi + n;
505 for(; miik < min_end;)
506 {
507 s33 += (*mikj) * (*(miik++));
508 mikj += n;
509 }
510 *(mij++) = s33;
511 }
512 for(j = 1; j <= ni; j++)
513 {
514 s34 = 0.0;
515 G4ErrorMatrixIter miik = mii + j;
516 G4ErrorMatrixIter mikij = mii + j * n + j;
517 for(G4int k = j; k <= ni; k++)
518 {
519 s34 += *mikij * (*(miik++));
520 mikij += n;
521 }
522 *(mii + j) = s34;
523 }
524 mi += n;
525 mii += (n + 1);
526 }
527 G4int nxch = ir[n];
528 if(nxch == 0)
529 return 0;
530 for(G4int mq = 1; mq <= nxch; mq++)
531 {
532 G4int k = nxch - mq + 1;
533 G4int ij = ir[k];
534 G4int i = ij >> 12;
535 G4int j = ij % 4096;
536 G4ErrorMatrixIter mki = m.begin() + i - 1;
537 G4ErrorMatrixIter mkj = m.begin() + j - 1;
538 for(k = 1; k <= n; k++)
539 {
540 // 2/24/05 David Sachs fix of improper swap bug that was present
541 // for many years:
542 G4double ti = *mki; // 2/24/05
543 *mki = *mkj;
544 *mkj = ti; // 2/24/05
545 mki += n;
546 mkj += n;
547 }
548 }
549 return 0;
550}
anonymous_namespace{G4QuasiElRatios.cc}::mi
const G4double mi
Definition: G4QuasiElRatios.cc:65

References error(), m, anonymous_namespace{G4QuasiElRatios.cc}::mi, CLHEP::detail::n, num_col(), and num_row().

Referenced by invert().

◆ error()

void G4ErrorMatrix::error ( const char *  s)
static

Definition at line 855 of file G4ErrorMatrix.cc.

856{
857 std::ostringstream message;
858 message << msg;
859 G4Exception("G4ErrorMatrix::error()", "GEANT4e-Error", FatalException,
860 message, "Exiting to System.");
861}
FatalException
@ FatalException
Definition: G4ExceptionSeverity.hh:61
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35

References FatalException, and G4Exception().

Referenced by determinant(), dfact_matrix(), dfinv_matrix(), G4ErrorMatrix(), G4ErrorSymMatrix::G4ErrorSymMatrix(), invert(), G4ErrorSymMatrix::sub(), and sub().

◆ inverse()

G4ErrorMatrix G4ErrorMatrix::inverse ( G4int &  ierr) const
inline

◆ invert()

void G4ErrorMatrix::invert ( G4int &  ierr)
virtual

Definition at line 684 of file G4ErrorMatrix.cc.

685{
686 if(ncol != nrow)
687 {
688 error("G4ErrorMatrix::invert: G4ErrorMatrix is not NxN");
689 }
690
691 static G4ThreadLocal G4int max_array = 20;
692 static G4ThreadLocal G4int* ir = 0;
693 if(!ir)
694 ir = new G4int[max_array + 1];
695
696 if(ncol > max_array)
697 {
698 delete[] ir;
699 max_array = nrow;
700 ir = new G4int[max_array + 1];
701 }
702 G4double t1, t2, t3;
703 G4double det, temp, ss;
704 G4int ifail;
705 switch(nrow)
706 {
707 case 3:
708 G4double c11, c12, c13, c21, c22, c23, c31, c32, c33;
709 ifail = 0;
710 c11 = (*(m.begin() + 4)) * (*(m.begin() + 8)) -
711 (*(m.begin() + 5)) * (*(m.begin() + 7));
712 c12 = (*(m.begin() + 5)) * (*(m.begin() + 6)) -
713 (*(m.begin() + 3)) * (*(m.begin() + 8));
714 c13 = (*(m.begin() + 3)) * (*(m.begin() + 7)) -
715 (*(m.begin() + 4)) * (*(m.begin() + 6));
716 c21 = (*(m.begin() + 7)) * (*(m.begin() + 2)) -
717 (*(m.begin() + 8)) * (*(m.begin() + 1));
718 c22 = (*(m.begin() + 8)) * (*m.begin()) -
719 (*(m.begin() + 6)) * (*(m.begin() + 2));
720 c23 = (*(m.begin() + 6)) * (*(m.begin() + 1)) -
721 (*(m.begin() + 7)) * (*m.begin());
722 c31 = (*(m.begin() + 1)) * (*(m.begin() + 5)) -
723 (*(m.begin() + 2)) * (*(m.begin() + 4));
724 c32 = (*(m.begin() + 2)) * (*(m.begin() + 3)) -
725 (*m.begin()) * (*(m.begin() + 5));
726 c33 = (*m.begin()) * (*(m.begin() + 4)) -
727 (*(m.begin() + 1)) * (*(m.begin() + 3));
728 t1 = std::fabs(*m.begin());
729 t2 = std::fabs(*(m.begin() + 3));
730 t3 = std::fabs(*(m.begin() + 6));
731 if(t1 >= t2)
732 {
733 if(t3 >= t1)
734 {
735 temp = *(m.begin() + 6);
736 det = c23 * c12 - c22 * c13;
737 }
738 else
739 {
740 temp = *(m.begin());
741 det = c22 * c33 - c23 * c32;
742 }
743 }
744 else if(t3 >= t2)
745 {
746 temp = *(m.begin() + 6);
747 det = c23 * c12 - c22 * c13;
748 }
749 else
750 {
751 temp = *(m.begin() + 3);
752 det = c13 * c32 - c12 * c33;
753 }
754 if(det == 0)
755 {
756 ierr = 1;
757 return;
758 }
759 {
760 ss = temp / det;
761 G4ErrorMatrixIter mq = m.begin();
762 *(mq++) = ss * c11;
763 *(mq++) = ss * c21;
764 *(mq++) = ss * c31;
765 *(mq++) = ss * c12;
766 *(mq++) = ss * c22;
767 *(mq++) = ss * c32;
768 *(mq++) = ss * c13;
769 *(mq++) = ss * c23;
770 *(mq) = ss * c33;
771 }
772 break;
773 case 2:
774 ifail = 0;
775 det = (*m.begin()) * (*(m.begin() + 3)) -
776 (*(m.begin() + 1)) * (*(m.begin() + 2));
777 if(det == 0)
778 {
779 ierr = 1;
780 return;
781 }
782 ss = 1.0 / det;
783 temp = ss * (*(m.begin() + 3));
784 *(m.begin() + 1) *= -ss;
785 *(m.begin() + 2) *= -ss;
786 *(m.begin() + 3) = ss * (*m.begin());
787 *(m.begin()) = temp;
788 break;
789 case 1:
790 ifail = 0;
791 if((*(m.begin())) == 0)
792 {
793 ierr = 1;
794 return;
795 }
796 *(m.begin()) = 1.0 / (*(m.begin()));
797 break;
798 case 4:
799 invertHaywood4(ierr);
800 return;
801 case 5:
802 invertHaywood5(ierr);
803 return;
804 case 6:
805 invertHaywood6(ierr);
806 return;
807 default:
808 ifail = dfact_matrix(det, ir);
809 if(ifail)
810 {
811 ierr = 1;
812 return;
813 }
814 dfinv_matrix(ir);
815 break;
816 }
817 ierr = 0;
818 return;
819}
G4ErrorMatrix::invertHaywood4
virtual void invertHaywood4(G4int &ierr)
Definition: G4ErrorMatrix.cc:959
G4ErrorMatrix::dfinv_matrix
G4int dfinv_matrix(G4int *ir)
Definition: G4ErrorMatrix.cc:432
G4ErrorMatrix::dfact_matrix
G4int dfact_matrix(G4double &det, G4int *ir)
Definition: G4ErrorMatrix.cc:552
G4ErrorMatrix::invertHaywood5
virtual void invertHaywood5(G4int &ierr)
Definition: G4ErrorMatrix.cc:1056
G4ErrorMatrix::invertHaywood6
virtual void invertHaywood6(G4int &ierr)
Definition: G4ErrorMatrix.cc:1277

References dfact_matrix(), dfinv_matrix(), error(), G4ThreadLocal, invertHaywood4(), invertHaywood5(), invertHaywood6(), m, ncol, and nrow.

◆ invertHaywood4()

void G4ErrorMatrix::invertHaywood4 ( G4int &  ierr)
protectedvirtual

Definition at line 959 of file G4ErrorMatrix.cc.

960{
961 ifail = 0;
962
963 // Find all NECESSARY 2x2 dets: (18 of them)
964
965 G4double Det2_12_01 = m[F10] * m[F21] - m[F11] * m[F20];
966 G4double Det2_12_02 = m[F10] * m[F22] - m[F12] * m[F20];
967 G4double Det2_12_03 = m[F10] * m[F23] - m[F13] * m[F20];
968 G4double Det2_12_13 = m[F11] * m[F23] - m[F13] * m[F21];
969 G4double Det2_12_23 = m[F12] * m[F23] - m[F13] * m[F22];
970 G4double Det2_12_12 = m[F11] * m[F22] - m[F12] * m[F21];
971 G4double Det2_13_01 = m[F10] * m[F31] - m[F11] * m[F30];
972 G4double Det2_13_02 = m[F10] * m[F32] - m[F12] * m[F30];
973 G4double Det2_13_03 = m[F10] * m[F33] - m[F13] * m[F30];
974 G4double Det2_13_12 = m[F11] * m[F32] - m[F12] * m[F31];
975 G4double Det2_13_13 = m[F11] * m[F33] - m[F13] * m[F31];
976 G4double Det2_13_23 = m[F12] * m[F33] - m[F13] * m[F32];
977 G4double Det2_23_01 = m[F20] * m[F31] - m[F21] * m[F30];
978 G4double Det2_23_02 = m[F20] * m[F32] - m[F22] * m[F30];
979 G4double Det2_23_03 = m[F20] * m[F33] - m[F23] * m[F30];
980 G4double Det2_23_12 = m[F21] * m[F32] - m[F22] * m[F31];
981 G4double Det2_23_13 = m[F21] * m[F33] - m[F23] * m[F31];
982 G4double Det2_23_23 = m[F22] * m[F33] - m[F23] * m[F32];
983
984 // Find all NECESSARY 3x3 dets: (16 of them)
985
986 G4double Det3_012_012 =
987 m[F00] * Det2_12_12 - m[F01] * Det2_12_02 + m[F02] * Det2_12_01;
988 G4double Det3_012_013 =
989 m[F00] * Det2_12_13 - m[F01] * Det2_12_03 + m[F03] * Det2_12_01;
990 G4double Det3_012_023 =
991 m[F00] * Det2_12_23 - m[F02] * Det2_12_03 + m[F03] * Det2_12_02;
992 G4double Det3_012_123 =
993 m[F01] * Det2_12_23 - m[F02] * Det2_12_13 + m[F03] * Det2_12_12;
994 G4double Det3_013_012 =
995 m[F00] * Det2_13_12 - m[F01] * Det2_13_02 + m[F02] * Det2_13_01;
996 G4double Det3_013_013 =
997 m[F00] * Det2_13_13 - m[F01] * Det2_13_03 + m[F03] * Det2_13_01;
998 G4double Det3_013_023 =
999 m[F00] * Det2_13_23 - m[F02] * Det2_13_03 + m[F03] * Det2_13_02;
1000 G4double Det3_013_123 =
1001 m[F01] * Det2_13_23 - m[F02] * Det2_13_13 + m[F03] * Det2_13_12;
1002 G4double Det3_023_012 =
1003 m[F00] * Det2_23_12 - m[F01] * Det2_23_02 + m[F02] * Det2_23_01;
1004 G4double Det3_023_013 =
1005 m[F00] * Det2_23_13 - m[F01] * Det2_23_03 + m[F03] * Det2_23_01;
1006 G4double Det3_023_023 =
1007 m[F00] * Det2_23_23 - m[F02] * Det2_23_03 + m[F03] * Det2_23_02;
1008 G4double Det3_023_123 =
1009 m[F01] * Det2_23_23 - m[F02] * Det2_23_13 + m[F03] * Det2_23_12;
1010 G4double Det3_123_012 =
1011 m[F10] * Det2_23_12 - m[F11] * Det2_23_02 + m[F12] * Det2_23_01;
1012 G4double Det3_123_013 =
1013 m[F10] * Det2_23_13 - m[F11] * Det2_23_03 + m[F13] * Det2_23_01;
1014 G4double Det3_123_023 =
1015 m[F10] * Det2_23_23 - m[F12] * Det2_23_03 + m[F13] * Det2_23_02;
1016 G4double Det3_123_123 =
1017 m[F11] * Det2_23_23 - m[F12] * Det2_23_13 + m[F13] * Det2_23_12;
1018
1019 // Find the 4x4 det:
1020
1021 G4double det = m[F00] * Det3_123_123 - m[F01] * Det3_123_023 +
1022 m[F02] * Det3_123_013 - m[F03] * Det3_123_012;
1023
1024 if(det == 0)
1025 {
1026 ifail = 1;
1027 return;
1028 }
1029
1030 G4double oneOverDet = 1.0 / det;
1031 G4double mn1OverDet = -oneOverDet;
1032
1033 m[F00] = Det3_123_123 * oneOverDet;
1034 m[F01] = Det3_023_123 * mn1OverDet;
1035 m[F02] = Det3_013_123 * oneOverDet;
1036 m[F03] = Det3_012_123 * mn1OverDet;
1037
1038 m[F10] = Det3_123_023 * mn1OverDet;
1039 m[F11] = Det3_023_023 * oneOverDet;
1040 m[F12] = Det3_013_023 * mn1OverDet;
1041 m[F13] = Det3_012_023 * oneOverDet;
1042
1043 m[F20] = Det3_123_013 * oneOverDet;
1044 m[F21] = Det3_023_013 * mn1OverDet;
1045 m[F22] = Det3_013_013 * oneOverDet;
1046 m[F23] = Det3_012_013 * mn1OverDet;
1047
1048 m[F30] = Det3_123_012 * mn1OverDet;
1049 m[F31] = Det3_023_012 * oneOverDet;
1050 m[F32] = Det3_013_012 * mn1OverDet;
1051 m[F33] = Det3_012_012 * oneOverDet;
1052
1053 return;
1054}
F01
#define F01
Definition: G4ErrorMatrix.cc:940
F22
#define F22
Definition: G4ErrorMatrix.cc:951
F00
#define F00
Definition: G4ErrorMatrix.cc:939
F02
#define F02
Definition: G4ErrorMatrix.cc:941
F31
#define F31
Definition: G4ErrorMatrix.cc:955
F32
#define F32
Definition: G4ErrorMatrix.cc:956
F21
#define F21
Definition: G4ErrorMatrix.cc:950
F03
#define F03
Definition: G4ErrorMatrix.cc:942
F30
#define F30
Definition: G4ErrorMatrix.cc:954
F11
#define F11
Definition: G4ErrorMatrix.cc:945
F10
#define F10
Definition: G4ErrorMatrix.cc:944
F33
#define F33
Definition: G4ErrorMatrix.cc:957
F20
#define F20
Definition: G4ErrorMatrix.cc:949
F12
#define F12
Definition: G4ErrorMatrix.cc:946
F23
#define F23
Definition: G4ErrorMatrix.cc:952
F13
#define F13
Definition: G4ErrorMatrix.cc:947

References F00, F01, F02, F03, F10, F11, F12, F13, F20, F21, F22, F23, F30, F31, F32, F33, and m.

Referenced by invert().

◆ invertHaywood5()

void G4ErrorMatrix::invertHaywood5 ( G4int &  ierr)
protectedvirtual

Definition at line 1056 of file G4ErrorMatrix.cc.

1057{
1058 ifail = 0;
1059
1060 // Find all NECESSARY 2x2 dets: (30 of them)
1061
1062 G4double Det2_23_01 = m[M20] * m[M31] - m[M21] * m[M30];
1063 G4double Det2_23_02 = m[M20] * m[M32] - m[M22] * m[M30];
1064 G4double Det2_23_03 = m[M20] * m[M33] - m[M23] * m[M30];
1065 G4double Det2_23_04 = m[M20] * m[M34] - m[M24] * m[M30];
1066 G4double Det2_23_12 = m[M21] * m[M32] - m[M22] * m[M31];
1067 G4double Det2_23_13 = m[M21] * m[M33] - m[M23] * m[M31];
1068 G4double Det2_23_14 = m[M21] * m[M34] - m[M24] * m[M31];
1069 G4double Det2_23_23 = m[M22] * m[M33] - m[M23] * m[M32];
1070 G4double Det2_23_24 = m[M22] * m[M34] - m[M24] * m[M32];
1071 G4double Det2_23_34 = m[M23] * m[M34] - m[M24] * m[M33];
1072 G4double Det2_24_01 = m[M20] * m[M41] - m[M21] * m[M40];
1073 G4double Det2_24_02 = m[M20] * m[M42] - m[M22] * m[M40];
1074 G4double Det2_24_03 = m[M20] * m[M43] - m[M23] * m[M40];
1075 G4double Det2_24_04 = m[M20] * m[M44] - m[M24] * m[M40];
1076 G4double Det2_24_12 = m[M21] * m[M42] - m[M22] * m[M41];
1077 G4double Det2_24_13 = m[M21] * m[M43] - m[M23] * m[M41];
1078 G4double Det2_24_14 = m[M21] * m[M44] - m[M24] * m[M41];
1079 G4double Det2_24_23 = m[M22] * m[M43] - m[M23] * m[M42];
1080 G4double Det2_24_24 = m[M22] * m[M44] - m[M24] * m[M42];
1081 G4double Det2_24_34 = m[M23] * m[M44] - m[M24] * m[M43];
1082 G4double Det2_34_01 = m[M30] * m[M41] - m[M31] * m[M40];
1083 G4double Det2_34_02 = m[M30] * m[M42] - m[M32] * m[M40];
1084 G4double Det2_34_03 = m[M30] * m[M43] - m[M33] * m[M40];
1085 G4double Det2_34_04 = m[M30] * m[M44] - m[M34] * m[M40];
1086 G4double Det2_34_12 = m[M31] * m[M42] - m[M32] * m[M41];
1087 G4double Det2_34_13 = m[M31] * m[M43] - m[M33] * m[M41];
1088 G4double Det2_34_14 = m[M31] * m[M44] - m[M34] * m[M41];
1089 G4double Det2_34_23 = m[M32] * m[M43] - m[M33] * m[M42];
1090 G4double Det2_34_24 = m[M32] * m[M44] - m[M34] * m[M42];
1091 G4double Det2_34_34 = m[M33] * m[M44] - m[M34] * m[M43];
1092
1093 // Find all NECESSARY 3x3 dets: (40 of them)
1094
1095 G4double Det3_123_012 =
1096 m[M10] * Det2_23_12 - m[M11] * Det2_23_02 + m[M12] * Det2_23_01;
1097 G4double Det3_123_013 =
1098 m[M10] * Det2_23_13 - m[M11] * Det2_23_03 + m[M13] * Det2_23_01;
1099 G4double Det3_123_014 =
1100 m[M10] * Det2_23_14 - m[M11] * Det2_23_04 + m[M14] * Det2_23_01;
1101 G4double Det3_123_023 =
1102 m[M10] * Det2_23_23 - m[M12] * Det2_23_03 + m[M13] * Det2_23_02;
1103 G4double Det3_123_024 =
1104 m[M10] * Det2_23_24 - m[M12] * Det2_23_04 + m[M14] * Det2_23_02;
1105 G4double Det3_123_034 =
1106 m[M10] * Det2_23_34 - m[M13] * Det2_23_04 + m[M14] * Det2_23_03;
1107 G4double Det3_123_123 =
1108 m[M11] * Det2_23_23 - m[M12] * Det2_23_13 + m[M13] * Det2_23_12;
1109 G4double Det3_123_124 =
1110 m[M11] * Det2_23_24 - m[M12] * Det2_23_14 + m[M14] * Det2_23_12;
1111 G4double Det3_123_134 =
1112 m[M11] * Det2_23_34 - m[M13] * Det2_23_14 + m[M14] * Det2_23_13;
1113 G4double Det3_123_234 =
1114 m[M12] * Det2_23_34 - m[M13] * Det2_23_24 + m[M14] * Det2_23_23;
1115 G4double Det3_124_012 =
1116 m[M10] * Det2_24_12 - m[M11] * Det2_24_02 + m[M12] * Det2_24_01;
1117 G4double Det3_124_013 =
1118 m[M10] * Det2_24_13 - m[M11] * Det2_24_03 + m[M13] * Det2_24_01;
1119 G4double Det3_124_014 =
1120 m[M10] * Det2_24_14 - m[M11] * Det2_24_04 + m[M14] * Det2_24_01;
1121 G4double Det3_124_023 =
1122 m[M10] * Det2_24_23 - m[M12] * Det2_24_03 + m[M13] * Det2_24_02;
1123 G4double Det3_124_024 =
1124 m[M10] * Det2_24_24 - m[M12] * Det2_24_04 + m[M14] * Det2_24_02;
1125 G4double Det3_124_034 =
1126 m[M10] * Det2_24_34 - m[M13] * Det2_24_04 + m[M14] * Det2_24_03;
1127 G4double Det3_124_123 =
1128 m[M11] * Det2_24_23 - m[M12] * Det2_24_13 + m[M13] * Det2_24_12;
1129 G4double Det3_124_124 =
1130 m[M11] * Det2_24_24 - m[M12] * Det2_24_14 + m[M14] * Det2_24_12;
1131 G4double Det3_124_134 =
1132 m[M11] * Det2_24_34 - m[M13] * Det2_24_14 + m[M14] * Det2_24_13;
1133 G4double Det3_124_234 =
1134 m[M12] * Det2_24_34 - m[M13] * Det2_24_24 + m[M14] * Det2_24_23;
1135 G4double Det3_134_012 =
1136 m[M10] * Det2_34_12 - m[M11] * Det2_34_02 + m[M12] * Det2_34_01;
1137 G4double Det3_134_013 =
1138 m[M10] * Det2_34_13 - m[M11] * Det2_34_03 + m[M13] * Det2_34_01;
1139 G4double Det3_134_014 =
1140 m[M10] * Det2_34_14 - m[M11] * Det2_34_04 + m[M14] * Det2_34_01;
1141 G4double Det3_134_023 =
1142 m[M10] * Det2_34_23 - m[M12] * Det2_34_03 + m[M13] * Det2_34_02;
1143 G4double Det3_134_024 =
1144 m[M10] * Det2_34_24 - m[M12] * Det2_34_04 + m[M14] * Det2_34_02;
1145 G4double Det3_134_034 =
1146 m[M10] * Det2_34_34 - m[M13] * Det2_34_04 + m[M14] * Det2_34_03;
1147 G4double Det3_134_123 =
1148 m[M11] * Det2_34_23 - m[M12] * Det2_34_13 + m[M13] * Det2_34_12;
1149 G4double Det3_134_124 =
1150 m[M11] * Det2_34_24 - m[M12] * Det2_34_14 + m[M14] * Det2_34_12;
1151 G4double Det3_134_134 =
1152 m[M11] * Det2_34_34 - m[M13] * Det2_34_14 + m[M14] * Det2_34_13;
1153 G4double Det3_134_234 =
1154 m[M12] * Det2_34_34 - m[M13] * Det2_34_24 + m[M14] * Det2_34_23;
1155 G4double Det3_234_012 =
1156 m[M20] * Det2_34_12 - m[M21] * Det2_34_02 + m[M22] * Det2_34_01;
1157 G4double Det3_234_013 =
1158 m[M20] * Det2_34_13 - m[M21] * Det2_34_03 + m[M23] * Det2_34_01;
1159 G4double Det3_234_014 =
1160 m[M20] * Det2_34_14 - m[M21] * Det2_34_04 + m[M24] * Det2_34_01;
1161 G4double Det3_234_023 =
1162 m[M20] * Det2_34_23 - m[M22] * Det2_34_03 + m[M23] * Det2_34_02;
1163 G4double Det3_234_024 =
1164 m[M20] * Det2_34_24 - m[M22] * Det2_34_04 + m[M24] * Det2_34_02;
1165 G4double Det3_234_034 =
1166 m[M20] * Det2_34_34 - m[M23] * Det2_34_04 + m[M24] * Det2_34_03;
1167 G4double Det3_234_123 =
1168 m[M21] * Det2_34_23 - m[M22] * Det2_34_13 + m[M23] * Det2_34_12;
1169 G4double Det3_234_124 =
1170 m[M21] * Det2_34_24 - m[M22] * Det2_34_14 + m[M24] * Det2_34_12;
1171 G4double Det3_234_134 =
1172 m[M21] * Det2_34_34 - m[M23] * Det2_34_14 + m[M24] * Det2_34_13;
1173 G4double Det3_234_234 =
1174 m[M22] * Det2_34_34 - m[M23] * Det2_34_24 + m[M24] * Det2_34_23;
1175
1176 // Find all NECESSARY 4x4 dets: (25 of them)
1177
1178 G4double Det4_0123_0123 = m[M00] * Det3_123_123 - m[M01] * Det3_123_023 +
1179 m[M02] * Det3_123_013 - m[M03] * Det3_123_012;
1180 G4double Det4_0123_0124 = m[M00] * Det3_123_124 - m[M01] * Det3_123_024 +
1181 m[M02] * Det3_123_014 - m[M04] * Det3_123_012;
1182 G4double Det4_0123_0134 = m[M00] * Det3_123_134 - m[M01] * Det3_123_034 +
1183 m[M03] * Det3_123_014 - m[M04] * Det3_123_013;
1184 G4double Det4_0123_0234 = m[M00] * Det3_123_234 - m[M02] * Det3_123_034 +
1185 m[M03] * Det3_123_024 - m[M04] * Det3_123_023;
1186 G4double Det4_0123_1234 = m[M01] * Det3_123_234 - m[M02] * Det3_123_134 +
1187 m[M03] * Det3_123_124 - m[M04] * Det3_123_123;
1188 G4double Det4_0124_0123 = m[M00] * Det3_124_123 - m[M01] * Det3_124_023 +
1189 m[M02] * Det3_124_013 - m[M03] * Det3_124_012;
1190 G4double Det4_0124_0124 = m[M00] * Det3_124_124 - m[M01] * Det3_124_024 +
1191 m[M02] * Det3_124_014 - m[M04] * Det3_124_012;
1192 G4double Det4_0124_0134 = m[M00] * Det3_124_134 - m[M01] * Det3_124_034 +
1193 m[M03] * Det3_124_014 - m[M04] * Det3_124_013;
1194 G4double Det4_0124_0234 = m[M00] * Det3_124_234 - m[M02] * Det3_124_034 +
1195 m[M03] * Det3_124_024 - m[M04] * Det3_124_023;
1196 G4double Det4_0124_1234 = m[M01] * Det3_124_234 - m[M02] * Det3_124_134 +
1197 m[M03] * Det3_124_124 - m[M04] * Det3_124_123;
1198 G4double Det4_0134_0123 = m[M00] * Det3_134_123 - m[M01] * Det3_134_023 +
1199 m[M02] * Det3_134_013 - m[M03] * Det3_134_012;
1200 G4double Det4_0134_0124 = m[M00] * Det3_134_124 - m[M01] * Det3_134_024 +
1201 m[M02] * Det3_134_014 - m[M04] * Det3_134_012;
1202 G4double Det4_0134_0134 = m[M00] * Det3_134_134 - m[M01] * Det3_134_034 +
1203 m[M03] * Det3_134_014 - m[M04] * Det3_134_013;
1204 G4double Det4_0134_0234 = m[M00] * Det3_134_234 - m[M02] * Det3_134_034 +
1205 m[M03] * Det3_134_024 - m[M04] * Det3_134_023;
1206 G4double Det4_0134_1234 = m[M01] * Det3_134_234 - m[M02] * Det3_134_134 +
1207 m[M03] * Det3_134_124 - m[M04] * Det3_134_123;
1208 G4double Det4_0234_0123 = m[M00] * Det3_234_123 - m[M01] * Det3_234_023 +
1209 m[M02] * Det3_234_013 - m[M03] * Det3_234_012;
1210 G4double Det4_0234_0124 = m[M00] * Det3_234_124 - m[M01] * Det3_234_024 +
1211 m[M02] * Det3_234_014 - m[M04] * Det3_234_012;
1212 G4double Det4_0234_0134 = m[M00] * Det3_234_134 - m[M01] * Det3_234_034 +
1213 m[M03] * Det3_234_014 - m[M04] * Det3_234_013;
1214 G4double Det4_0234_0234 = m[M00] * Det3_234_234 - m[M02] * Det3_234_034 +
1215 m[M03] * Det3_234_024 - m[M04] * Det3_234_023;
1216 G4double Det4_0234_1234 = m[M01] * Det3_234_234 - m[M02] * Det3_234_134 +
1217 m[M03] * Det3_234_124 - m[M04] * Det3_234_123;
1218 G4double Det4_1234_0123 = m[M10] * Det3_234_123 - m[M11] * Det3_234_023 +
1219 m[M12] * Det3_234_013 - m[M13] * Det3_234_012;
1220 G4double Det4_1234_0124 = m[M10] * Det3_234_124 - m[M11] * Det3_234_024 +
1221 m[M12] * Det3_234_014 - m[M14] * Det3_234_012;
1222 G4double Det4_1234_0134 = m[M10] * Det3_234_134 - m[M11] * Det3_234_034 +
1223 m[M13] * Det3_234_014 - m[M14] * Det3_234_013;
1224 G4double Det4_1234_0234 = m[M10] * Det3_234_234 - m[M12] * Det3_234_034 +
1225 m[M13] * Det3_234_024 - m[M14] * Det3_234_023;
1226 G4double Det4_1234_1234 = m[M11] * Det3_234_234 - m[M12] * Det3_234_134 +
1227 m[M13] * Det3_234_124 - m[M14] * Det3_234_123;
1228
1229 // Find the 5x5 det:
1230
1231 G4double det = m[M00] * Det4_1234_1234 - m[M01] * Det4_1234_0234 +
1232 m[M02] * Det4_1234_0134 - m[M03] * Det4_1234_0124 +
1233 m[M04] * Det4_1234_0123;
1234
1235 if(det == 0)
1236 {
1237 ifail = 1;
1238 return;
1239 }
1240
1241 G4double oneOverDet = 1.0 / det;
1242 G4double mn1OverDet = -oneOverDet;
1243
1244 m[M00] = Det4_1234_1234 * oneOverDet;
1245 m[M01] = Det4_0234_1234 * mn1OverDet;
1246 m[M02] = Det4_0134_1234 * oneOverDet;
1247 m[M03] = Det4_0124_1234 * mn1OverDet;
1248 m[M04] = Det4_0123_1234 * oneOverDet;
1249
1250 m[M10] = Det4_1234_0234 * mn1OverDet;
1251 m[M11] = Det4_0234_0234 * oneOverDet;
1252 m[M12] = Det4_0134_0234 * mn1OverDet;
1253 m[M13] = Det4_0124_0234 * oneOverDet;
1254 m[M14] = Det4_0123_0234 * mn1OverDet;
1255
1256 m[M20] = Det4_1234_0134 * oneOverDet;
1257 m[M21] = Det4_0234_0134 * mn1OverDet;
1258 m[M22] = Det4_0134_0134 * oneOverDet;
1259 m[M23] = Det4_0124_0134 * mn1OverDet;
1260 m[M24] = Det4_0123_0134 * oneOverDet;
1261
1262 m[M30] = Det4_1234_0124 * mn1OverDet;
1263 m[M31] = Det4_0234_0124 * oneOverDet;
1264 m[M32] = Det4_0134_0124 * mn1OverDet;
1265 m[M33] = Det4_0124_0124 * oneOverDet;
1266 m[M34] = Det4_0123_0124 * mn1OverDet;
1267
1268 m[M40] = Det4_1234_0123 * oneOverDet;
1269 m[M41] = Det4_0234_0123 * mn1OverDet;
1270 m[M42] = Det4_0134_0123 * oneOverDet;
1271 m[M43] = Det4_0124_0123 * mn1OverDet;
1272 m[M44] = Det4_0123_0123 * oneOverDet;
1273
1274 return;
1275}
M00
#define M00
Definition: G4ErrorMatrix.cc:909
M34
#define M34
Definition: G4ErrorMatrix.cc:931
M03
#define M03
Definition: G4ErrorMatrix.cc:912
M10
#define M10
Definition: G4ErrorMatrix.cc:915
M41
#define M41
Definition: G4ErrorMatrix.cc:934
M20
#define M20
Definition: G4ErrorMatrix.cc:921
M42
#define M42
Definition: G4ErrorMatrix.cc:935
M13
#define M13
Definition: G4ErrorMatrix.cc:918
M33
#define M33
Definition: G4ErrorMatrix.cc:930
M04
#define M04
Definition: G4ErrorMatrix.cc:913
M23
#define M23
Definition: G4ErrorMatrix.cc:924
M11
#define M11
Definition: G4ErrorMatrix.cc:916
M44
#define M44
Definition: G4ErrorMatrix.cc:937
M31
#define M31
Definition: G4ErrorMatrix.cc:928
M02
#define M02
Definition: G4ErrorMatrix.cc:911
M21
#define M21
Definition: G4ErrorMatrix.cc:922
M01
#define M01
Definition: G4ErrorMatrix.cc:910
M24
#define M24
Definition: G4ErrorMatrix.cc:925
M12
#define M12
Definition: G4ErrorMatrix.cc:917
M22
#define M22
Definition: G4ErrorMatrix.cc:923
M43
#define M43
Definition: G4ErrorMatrix.cc:936
M14
#define M14
Definition: G4ErrorMatrix.cc:919
M30
#define M30
Definition: G4ErrorMatrix.cc:927
M40
#define M40
Definition: G4ErrorMatrix.cc:933
M32
#define M32
Definition: G4ErrorMatrix.cc:929

References m, M00, M01, M02, M03, M04, M10, M11, M12, M13, M14, M20, M21, M22, M23, M24, M30, M31, M32, M33, M34, M40, M41, M42, M43, and M44.

Referenced by invert().

◆ invertHaywood6()

void G4ErrorMatrix::invertHaywood6 ( G4int &  ierr)
protectedvirtual

Definition at line 1277 of file G4ErrorMatrix.cc.

1278{
1279 ifail = 0;
1280
1281 // Find all NECESSARY 2x2 dets: (45 of them)
1282
1283 G4double Det2_34_01 = m[A30] * m[A41] - m[A31] * m[A40];
1284 G4double Det2_34_02 = m[A30] * m[A42] - m[A32] * m[A40];
1285 G4double Det2_34_03 = m[A30] * m[A43] - m[A33] * m[A40];
1286 G4double Det2_34_04 = m[A30] * m[A44] - m[A34] * m[A40];
1287 G4double Det2_34_05 = m[A30] * m[A45] - m[A35] * m[A40];
1288 G4double Det2_34_12 = m[A31] * m[A42] - m[A32] * m[A41];
1289 G4double Det2_34_13 = m[A31] * m[A43] - m[A33] * m[A41];
1290 G4double Det2_34_14 = m[A31] * m[A44] - m[A34] * m[A41];
1291 G4double Det2_34_15 = m[A31] * m[A45] - m[A35] * m[A41];
1292 G4double Det2_34_23 = m[A32] * m[A43] - m[A33] * m[A42];
1293 G4double Det2_34_24 = m[A32] * m[A44] - m[A34] * m[A42];
1294 G4double Det2_34_25 = m[A32] * m[A45] - m[A35] * m[A42];
1295 G4double Det2_34_34 = m[A33] * m[A44] - m[A34] * m[A43];
1296 G4double Det2_34_35 = m[A33] * m[A45] - m[A35] * m[A43];
1297 G4double Det2_34_45 = m[A34] * m[A45] - m[A35] * m[A44];
1298 G4double Det2_35_01 = m[A30] * m[A51] - m[A31] * m[A50];
1299 G4double Det2_35_02 = m[A30] * m[A52] - m[A32] * m[A50];
1300 G4double Det2_35_03 = m[A30] * m[A53] - m[A33] * m[A50];
1301 G4double Det2_35_04 = m[A30] * m[A54] - m[A34] * m[A50];
1302 G4double Det2_35_05 = m[A30] * m[A55] - m[A35] * m[A50];
1303 G4double Det2_35_12 = m[A31] * m[A52] - m[A32] * m[A51];
1304 G4double Det2_35_13 = m[A31] * m[A53] - m[A33] * m[A51];
1305 G4double Det2_35_14 = m[A31] * m[A54] - m[A34] * m[A51];
1306 G4double Det2_35_15 = m[A31] * m[A55] - m[A35] * m[A51];
1307 G4double Det2_35_23 = m[A32] * m[A53] - m[A33] * m[A52];
1308 G4double Det2_35_24 = m[A32] * m[A54] - m[A34] * m[A52];
1309 G4double Det2_35_25 = m[A32] * m[A55] - m[A35] * m[A52];
1310 G4double Det2_35_34 = m[A33] * m[A54] - m[A34] * m[A53];
1311 G4double Det2_35_35 = m[A33] * m[A55] - m[A35] * m[A53];
1312 G4double Det2_35_45 = m[A34] * m[A55] - m[A35] * m[A54];
1313 G4double Det2_45_01 = m[A40] * m[A51] - m[A41] * m[A50];
1314 G4double Det2_45_02 = m[A40] * m[A52] - m[A42] * m[A50];
1315 G4double Det2_45_03 = m[A40] * m[A53] - m[A43] * m[A50];
1316 G4double Det2_45_04 = m[A40] * m[A54] - m[A44] * m[A50];
1317 G4double Det2_45_05 = m[A40] * m[A55] - m[A45] * m[A50];
1318 G4double Det2_45_12 = m[A41] * m[A52] - m[A42] * m[A51];
1319 G4double Det2_45_13 = m[A41] * m[A53] - m[A43] * m[A51];
1320 G4double Det2_45_14 = m[A41] * m[A54] - m[A44] * m[A51];
1321 G4double Det2_45_15 = m[A41] * m[A55] - m[A45] * m[A51];
1322 G4double Det2_45_23 = m[A42] * m[A53] - m[A43] * m[A52];
1323 G4double Det2_45_24 = m[A42] * m[A54] - m[A44] * m[A52];
1324 G4double Det2_45_25 = m[A42] * m[A55] - m[A45] * m[A52];
1325 G4double Det2_45_34 = m[A43] * m[A54] - m[A44] * m[A53];
1326 G4double Det2_45_35 = m[A43] * m[A55] - m[A45] * m[A53];
1327 G4double Det2_45_45 = m[A44] * m[A55] - m[A45] * m[A54];
1328
1329 // Find all NECESSARY 3x3 dets: (80 of them)
1330
1331 G4double Det3_234_012 =
1332 m[A20] * Det2_34_12 - m[A21] * Det2_34_02 + m[A22] * Det2_34_01;
1333 G4double Det3_234_013 =
1334 m[A20] * Det2_34_13 - m[A21] * Det2_34_03 + m[A23] * Det2_34_01;
1335 G4double Det3_234_014 =
1336 m[A20] * Det2_34_14 - m[A21] * Det2_34_04 + m[A24] * Det2_34_01;
1337 G4double Det3_234_015 =
1338 m[A20] * Det2_34_15 - m[A21] * Det2_34_05 + m[A25] * Det2_34_01;
1339 G4double Det3_234_023 =
1340 m[A20] * Det2_34_23 - m[A22] * Det2_34_03 + m[A23] * Det2_34_02;
1341 G4double Det3_234_024 =
1342 m[A20] * Det2_34_24 - m[A22] * Det2_34_04 + m[A24] * Det2_34_02;
1343 G4double Det3_234_025 =
1344 m[A20] * Det2_34_25 - m[A22] * Det2_34_05 + m[A25] * Det2_34_02;
1345 G4double Det3_234_034 =
1346 m[A20] * Det2_34_34 - m[A23] * Det2_34_04 + m[A24] * Det2_34_03;
1347 G4double Det3_234_035 =
1348 m[A20] * Det2_34_35 - m[A23] * Det2_34_05 + m[A25] * Det2_34_03;
1349 G4double Det3_234_045 =
1350 m[A20] * Det2_34_45 - m[A24] * Det2_34_05 + m[A25] * Det2_34_04;
1351 G4double Det3_234_123 =
1352 m[A21] * Det2_34_23 - m[A22] * Det2_34_13 + m[A23] * Det2_34_12;
1353 G4double Det3_234_124 =
1354 m[A21] * Det2_34_24 - m[A22] * Det2_34_14 + m[A24] * Det2_34_12;
1355 G4double Det3_234_125 =
1356 m[A21] * Det2_34_25 - m[A22] * Det2_34_15 + m[A25] * Det2_34_12;
1357 G4double Det3_234_134 =
1358 m[A21] * Det2_34_34 - m[A23] * Det2_34_14 + m[A24] * Det2_34_13;
1359 G4double Det3_234_135 =
1360 m[A21] * Det2_34_35 - m[A23] * Det2_34_15 + m[A25] * Det2_34_13;
1361 G4double Det3_234_145 =
1362 m[A21] * Det2_34_45 - m[A24] * Det2_34_15 + m[A25] * Det2_34_14;
1363 G4double Det3_234_234 =
1364 m[A22] * Det2_34_34 - m[A23] * Det2_34_24 + m[A24] * Det2_34_23;
1365 G4double Det3_234_235 =
1366 m[A22] * Det2_34_35 - m[A23] * Det2_34_25 + m[A25] * Det2_34_23;
1367 G4double Det3_234_245 =
1368 m[A22] * Det2_34_45 - m[A24] * Det2_34_25 + m[A25] * Det2_34_24;
1369 G4double Det3_234_345 =
1370 m[A23] * Det2_34_45 - m[A24] * Det2_34_35 + m[A25] * Det2_34_34;
1371 G4double Det3_235_012 =
1372 m[A20] * Det2_35_12 - m[A21] * Det2_35_02 + m[A22] * Det2_35_01;
1373 G4double Det3_235_013 =
1374 m[A20] * Det2_35_13 - m[A21] * Det2_35_03 + m[A23] * Det2_35_01;
1375 G4double Det3_235_014 =
1376 m[A20] * Det2_35_14 - m[A21] * Det2_35_04 + m[A24] * Det2_35_01;
1377 G4double Det3_235_015 =
1378 m[A20] * Det2_35_15 - m[A21] * Det2_35_05 + m[A25] * Det2_35_01;
1379 G4double Det3_235_023 =
1380 m[A20] * Det2_35_23 - m[A22] * Det2_35_03 + m[A23] * Det2_35_02;
1381 G4double Det3_235_024 =
1382 m[A20] * Det2_35_24 - m[A22] * Det2_35_04 + m[A24] * Det2_35_02;
1383 G4double Det3_235_025 =
1384 m[A20] * Det2_35_25 - m[A22] * Det2_35_05 + m[A25] * Det2_35_02;
1385 G4double Det3_235_034 =
1386 m[A20] * Det2_35_34 - m[A23] * Det2_35_04 + m[A24] * Det2_35_03;
1387 G4double Det3_235_035 =
1388 m[A20] * Det2_35_35 - m[A23] * Det2_35_05 + m[A25] * Det2_35_03;
1389 G4double Det3_235_045 =
1390 m[A20] * Det2_35_45 - m[A24] * Det2_35_05 + m[A25] * Det2_35_04;
1391 G4double Det3_235_123 =
1392 m[A21] * Det2_35_23 - m[A22] * Det2_35_13 + m[A23] * Det2_35_12;
1393 G4double Det3_235_124 =
1394 m[A21] * Det2_35_24 - m[A22] * Det2_35_14 + m[A24] * Det2_35_12;
1395 G4double Det3_235_125 =
1396 m[A21] * Det2_35_25 - m[A22] * Det2_35_15 + m[A25] * Det2_35_12;
1397 G4double Det3_235_134 =
1398 m[A21] * Det2_35_34 - m[A23] * Det2_35_14 + m[A24] * Det2_35_13;
1399 G4double Det3_235_135 =
1400 m[A21] * Det2_35_35 - m[A23] * Det2_35_15 + m[A25] * Det2_35_13;
1401 G4double Det3_235_145 =
1402 m[A21] * Det2_35_45 - m[A24] * Det2_35_15 + m[A25] * Det2_35_14;
1403 G4double Det3_235_234 =
1404 m[A22] * Det2_35_34 - m[A23] * Det2_35_24 + m[A24] * Det2_35_23;
1405 G4double Det3_235_235 =
1406 m[A22] * Det2_35_35 - m[A23] * Det2_35_25 + m[A25] * Det2_35_23;
1407 G4double Det3_235_245 =
1408 m[A22] * Det2_35_45 - m[A24] * Det2_35_25 + m[A25] * Det2_35_24;
1409 G4double Det3_235_345 =
1410 m[A23] * Det2_35_45 - m[A24] * Det2_35_35 + m[A25] * Det2_35_34;
1411 G4double Det3_245_012 =
1412 m[A20] * Det2_45_12 - m[A21] * Det2_45_02 + m[A22] * Det2_45_01;
1413 G4double Det3_245_013 =
1414 m[A20] * Det2_45_13 - m[A21] * Det2_45_03 + m[A23] * Det2_45_01;
1415 G4double Det3_245_014 =
1416 m[A20] * Det2_45_14 - m[A21] * Det2_45_04 + m[A24] * Det2_45_01;
1417 G4double Det3_245_015 =
1418 m[A20] * Det2_45_15 - m[A21] * Det2_45_05 + m[A25] * Det2_45_01;
1419 G4double Det3_245_023 =
1420 m[A20] * Det2_45_23 - m[A22] * Det2_45_03 + m[A23] * Det2_45_02;
1421 G4double Det3_245_024 =
1422 m[A20] * Det2_45_24 - m[A22] * Det2_45_04 + m[A24] * Det2_45_02;
1423 G4double Det3_245_025 =
1424 m[A20] * Det2_45_25 - m[A22] * Det2_45_05 + m[A25] * Det2_45_02;
1425 G4double Det3_245_034 =
1426 m[A20] * Det2_45_34 - m[A23] * Det2_45_04 + m[A24] * Det2_45_03;
1427 G4double Det3_245_035 =
1428 m[A20] * Det2_45_35 - m[A23] * Det2_45_05 + m[A25] * Det2_45_03;
1429 G4double Det3_245_045 =
1430 m[A20] * Det2_45_45 - m[A24] * Det2_45_05 + m[A25] * Det2_45_04;
1431 G4double Det3_245_123 =
1432 m[A21] * Det2_45_23 - m[A22] * Det2_45_13 + m[A23] * Det2_45_12;
1433 G4double Det3_245_124 =
1434 m[A21] * Det2_45_24 - m[A22] * Det2_45_14 + m[A24] * Det2_45_12;
1435 G4double Det3_245_125 =
1436 m[A21] * Det2_45_25 - m[A22] * Det2_45_15 + m[A25] * Det2_45_12;
1437 G4double Det3_245_134 =
1438 m[A21] * Det2_45_34 - m[A23] * Det2_45_14 + m[A24] * Det2_45_13;
1439 G4double Det3_245_135 =
1440 m[A21] * Det2_45_35 - m[A23] * Det2_45_15 + m[A25] * Det2_45_13;
1441 G4double Det3_245_145 =
1442 m[A21] * Det2_45_45 - m[A24] * Det2_45_15 + m[A25] * Det2_45_14;
1443 G4double Det3_245_234 =
1444 m[A22] * Det2_45_34 - m[A23] * Det2_45_24 + m[A24] * Det2_45_23;
1445 G4double Det3_245_235 =
1446 m[A22] * Det2_45_35 - m[A23] * Det2_45_25 + m[A25] * Det2_45_23;
1447 G4double Det3_245_245 =
1448 m[A22] * Det2_45_45 - m[A24] * Det2_45_25 + m[A25] * Det2_45_24;
1449 G4double Det3_245_345 =
1450 m[A23] * Det2_45_45 - m[A24] * Det2_45_35 + m[A25] * Det2_45_34;
1451 G4double Det3_345_012 =
1452 m[A30] * Det2_45_12 - m[A31] * Det2_45_02 + m[A32] * Det2_45_01;
1453 G4double Det3_345_013 =
1454 m[A30] * Det2_45_13 - m[A31] * Det2_45_03 + m[A33] * Det2_45_01;
1455 G4double Det3_345_014 =
1456 m[A30] * Det2_45_14 - m[A31] * Det2_45_04 + m[A34] * Det2_45_01;
1457 G4double Det3_345_015 =
1458 m[A30] * Det2_45_15 - m[A31] * Det2_45_05 + m[A35] * Det2_45_01;
1459 G4double Det3_345_023 =
1460 m[A30] * Det2_45_23 - m[A32] * Det2_45_03 + m[A33] * Det2_45_02;
1461 G4double Det3_345_024 =
1462 m[A30] * Det2_45_24 - m[A32] * Det2_45_04 + m[A34] * Det2_45_02;
1463 G4double Det3_345_025 =
1464 m[A30] * Det2_45_25 - m[A32] * Det2_45_05 + m[A35] * Det2_45_02;
1465 G4double Det3_345_034 =
1466 m[A30] * Det2_45_34 - m[A33] * Det2_45_04 + m[A34] * Det2_45_03;
1467 G4double Det3_345_035 =
1468 m[A30] * Det2_45_35 - m[A33] * Det2_45_05 + m[A35] * Det2_45_03;
1469 G4double Det3_345_045 =
1470 m[A30] * Det2_45_45 - m[A34] * Det2_45_05 + m[A35] * Det2_45_04;
1471 G4double Det3_345_123 =
1472 m[A31] * Det2_45_23 - m[A32] * Det2_45_13 + m[A33] * Det2_45_12;
1473 G4double Det3_345_124 =
1474 m[A31] * Det2_45_24 - m[A32] * Det2_45_14 + m[A34] * Det2_45_12;
1475 G4double Det3_345_125 =
1476 m[A31] * Det2_45_25 - m[A32] * Det2_45_15 + m[A35] * Det2_45_12;
1477 G4double Det3_345_134 =
1478 m[A31] * Det2_45_34 - m[A33] * Det2_45_14 + m[A34] * Det2_45_13;
1479 G4double Det3_345_135 =
1480 m[A31] * Det2_45_35 - m[A33] * Det2_45_15 + m[A35] * Det2_45_13;
1481 G4double Det3_345_145 =
1482 m[A31] * Det2_45_45 - m[A34] * Det2_45_15 + m[A35] * Det2_45_14;
1483 G4double Det3_345_234 =
1484 m[A32] * Det2_45_34 - m[A33] * Det2_45_24 + m[A34] * Det2_45_23;
1485 G4double Det3_345_235 =
1486 m[A32] * Det2_45_35 - m[A33] * Det2_45_25 + m[A35] * Det2_45_23;
1487 G4double Det3_345_245 =
1488 m[A32] * Det2_45_45 - m[A34] * Det2_45_25 + m[A35] * Det2_45_24;
1489 G4double Det3_345_345 =
1490 m[A33] * Det2_45_45 - m[A34] * Det2_45_35 + m[A35] * Det2_45_34;
1491
1492 // Find all NECESSARY 4x4 dets: (75 of them)
1493
1494 G4double Det4_1234_0123 = m[A10] * Det3_234_123 - m[A11] * Det3_234_023 +
1495 m[A12] * Det3_234_013 - m[A13] * Det3_234_012;
1496 G4double Det4_1234_0124 = m[A10] * Det3_234_124 - m[A11] * Det3_234_024 +
1497 m[A12] * Det3_234_014 - m[A14] * Det3_234_012;
1498 G4double Det4_1234_0125 = m[A10] * Det3_234_125 - m[A11] * Det3_234_025 +
1499 m[A12] * Det3_234_015 - m[A15] * Det3_234_012;
1500 G4double Det4_1234_0134 = m[A10] * Det3_234_134 - m[A11] * Det3_234_034 +
1501 m[A13] * Det3_234_014 - m[A14] * Det3_234_013;
1502 G4double Det4_1234_0135 = m[A10] * Det3_234_135 - m[A11] * Det3_234_035 +
1503 m[A13] * Det3_234_015 - m[A15] * Det3_234_013;
1504 G4double Det4_1234_0145 = m[A10] * Det3_234_145 - m[A11] * Det3_234_045 +
1505 m[A14] * Det3_234_015 - m[A15] * Det3_234_014;
1506 G4double Det4_1234_0234 = m[A10] * Det3_234_234 - m[A12] * Det3_234_034 +
1507 m[A13] * Det3_234_024 - m[A14] * Det3_234_023;
1508 G4double Det4_1234_0235 = m[A10] * Det3_234_235 - m[A12] * Det3_234_035 +
1509 m[A13] * Det3_234_025 - m[A15] * Det3_234_023;
1510 G4double Det4_1234_0245 = m[A10] * Det3_234_245 - m[A12] * Det3_234_045 +
1511 m[A14] * Det3_234_025 - m[A15] * Det3_234_024;
1512 G4double Det4_1234_0345 = m[A10] * Det3_234_345 - m[A13] * Det3_234_045 +
1513 m[A14] * Det3_234_035 - m[A15] * Det3_234_034;
1514 G4double Det4_1234_1234 = m[A11] * Det3_234_234 - m[A12] * Det3_234_134 +
1515 m[A13] * Det3_234_124 - m[A14] * Det3_234_123;
1516 G4double Det4_1234_1235 = m[A11] * Det3_234_235 - m[A12] * Det3_234_135 +
1517 m[A13] * Det3_234_125 - m[A15] * Det3_234_123;
1518 G4double Det4_1234_1245 = m[A11] * Det3_234_245 - m[A12] * Det3_234_145 +
1519 m[A14] * Det3_234_125 - m[A15] * Det3_234_124;
1520 G4double Det4_1234_1345 = m[A11] * Det3_234_345 - m[A13] * Det3_234_145 +
1521 m[A14] * Det3_234_135 - m[A15] * Det3_234_134;
1522 G4double Det4_1234_2345 = m[A12] * Det3_234_345 - m[A13] * Det3_234_245 +
1523 m[A14] * Det3_234_235 - m[A15] * Det3_234_234;
1524 G4double Det4_1235_0123 = m[A10] * Det3_235_123 - m[A11] * Det3_235_023 +
1525 m[A12] * Det3_235_013 - m[A13] * Det3_235_012;
1526 G4double Det4_1235_0124 = m[A10] * Det3_235_124 - m[A11] * Det3_235_024 +
1527 m[A12] * Det3_235_014 - m[A14] * Det3_235_012;
1528 G4double Det4_1235_0125 = m[A10] * Det3_235_125 - m[A11] * Det3_235_025 +
1529 m[A12] * Det3_235_015 - m[A15] * Det3_235_012;
1530 G4double Det4_1235_0134 = m[A10] * Det3_235_134 - m[A11] * Det3_235_034 +
1531 m[A13] * Det3_235_014 - m[A14] * Det3_235_013;
1532 G4double Det4_1235_0135 = m[A10] * Det3_235_135 - m[A11] * Det3_235_035 +
1533 m[A13] * Det3_235_015 - m[A15] * Det3_235_013;
1534 G4double Det4_1235_0145 = m[A10] * Det3_235_145 - m[A11] * Det3_235_045 +
1535 m[A14] * Det3_235_015 - m[A15] * Det3_235_014;
1536 G4double Det4_1235_0234 = m[A10] * Det3_235_234 - m[A12] * Det3_235_034 +
1537 m[A13] * Det3_235_024 - m[A14] * Det3_235_023;
1538 G4double Det4_1235_0235 = m[A10] * Det3_235_235 - m[A12] * Det3_235_035 +
1539 m[A13] * Det3_235_025 - m[A15] * Det3_235_023;
1540 G4double Det4_1235_0245 = m[A10] * Det3_235_245 - m[A12] * Det3_235_045 +
1541 m[A14] * Det3_235_025 - m[A15] * Det3_235_024;
1542 G4double Det4_1235_0345 = m[A10] * Det3_235_345 - m[A13] * Det3_235_045 +
1543 m[A14] * Det3_235_035 - m[A15] * Det3_235_034;
1544 G4double Det4_1235_1234 = m[A11] * Det3_235_234 - m[A12] * Det3_235_134 +
1545 m[A13] * Det3_235_124 - m[A14] * Det3_235_123;
1546 G4double Det4_1235_1235 = m[A11] * Det3_235_235 - m[A12] * Det3_235_135 +
1547 m[A13] * Det3_235_125 - m[A15] * Det3_235_123;
1548 G4double Det4_1235_1245 = m[A11] * Det3_235_245 - m[A12] * Det3_235_145 +
1549 m[A14] * Det3_235_125 - m[A15] * Det3_235_124;
1550 G4double Det4_1235_1345 = m[A11] * Det3_235_345 - m[A13] * Det3_235_145 +
1551 m[A14] * Det3_235_135 - m[A15] * Det3_235_134;
1552 G4double Det4_1235_2345 = m[A12] * Det3_235_345 - m[A13] * Det3_235_245 +
1553 m[A14] * Det3_235_235 - m[A15] * Det3_235_234;
1554 G4double Det4_1245_0123 = m[A10] * Det3_245_123 - m[A11] * Det3_245_023 +
1555 m[A12] * Det3_245_013 - m[A13] * Det3_245_012;
1556 G4double Det4_1245_0124 = m[A10] * Det3_245_124 - m[A11] * Det3_245_024 +
1557 m[A12] * Det3_245_014 - m[A14] * Det3_245_012;
1558 G4double Det4_1245_0125 = m[A10] * Det3_245_125 - m[A11] * Det3_245_025 +
1559 m[A12] * Det3_245_015 - m[A15] * Det3_245_012;
1560 G4double Det4_1245_0134 = m[A10] * Det3_245_134 - m[A11] * Det3_245_034 +
1561 m[A13] * Det3_245_014 - m[A14] * Det3_245_013;
1562 G4double Det4_1245_0135 = m[A10] * Det3_245_135 - m[A11] * Det3_245_035 +
1563 m[A13] * Det3_245_015 - m[A15] * Det3_245_013;
1564 G4double Det4_1245_0145 = m[A10] * Det3_245_145 - m[A11] * Det3_245_045 +
1565 m[A14] * Det3_245_015 - m[A15] * Det3_245_014;
1566 G4double Det4_1245_0234 = m[A10] * Det3_245_234 - m[A12] * Det3_245_034 +
1567 m[A13] * Det3_245_024 - m[A14] * Det3_245_023;
1568 G4double Det4_1245_0235 = m[A10] * Det3_245_235 - m[A12] * Det3_245_035 +
1569 m[A13] * Det3_245_025 - m[A15] * Det3_245_023;
1570 G4double Det4_1245_0245 = m[A10] * Det3_245_245 - m[A12] * Det3_245_045 +
1571 m[A14] * Det3_245_025 - m[A15] * Det3_245_024;
1572 G4double Det4_1245_0345 = m[A10] * Det3_245_345 - m[A13] * Det3_245_045 +
1573 m[A14] * Det3_245_035 - m[A15] * Det3_245_034;
1574 G4double Det4_1245_1234 = m[A11] * Det3_245_234 - m[A12] * Det3_245_134 +
1575 m[A13] * Det3_245_124 - m[A14] * Det3_245_123;
1576 G4double Det4_1245_1235 = m[A11] * Det3_245_235 - m[A12] * Det3_245_135 +
1577 m[A13] * Det3_245_125 - m[A15] * Det3_245_123;
1578 G4double Det4_1245_1245 = m[A11] * Det3_245_245 - m[A12] * Det3_245_145 +
1579 m[A14] * Det3_245_125 - m[A15] * Det3_245_124;
1580 G4double Det4_1245_1345 = m[A11] * Det3_245_345 - m[A13] * Det3_245_145 +
1581 m[A14] * Det3_245_135 - m[A15] * Det3_245_134;
1582 G4double Det4_1245_2345 = m[A12] * Det3_245_345 - m[A13] * Det3_245_245 +
1583 m[A14] * Det3_245_235 - m[A15] * Det3_245_234;
1584 G4double Det4_1345_0123 = m[A10] * Det3_345_123 - m[A11] * Det3_345_023 +
1585 m[A12] * Det3_345_013 - m[A13] * Det3_345_012;
1586 G4double Det4_1345_0124 = m[A10] * Det3_345_124 - m[A11] * Det3_345_024 +
1587 m[A12] * Det3_345_014 - m[A14] * Det3_345_012;
1588 G4double Det4_1345_0125 = m[A10] * Det3_345_125 - m[A11] * Det3_345_025 +
1589 m[A12] * Det3_345_015 - m[A15] * Det3_345_012;
1590 G4double Det4_1345_0134 = m[A10] * Det3_345_134 - m[A11] * Det3_345_034 +
1591 m[A13] * Det3_345_014 - m[A14] * Det3_345_013;
1592 G4double Det4_1345_0135 = m[A10] * Det3_345_135 - m[A11] * Det3_345_035 +
1593 m[A13] * Det3_345_015 - m[A15] * Det3_345_013;
1594 G4double Det4_1345_0145 = m[A10] * Det3_345_145 - m[A11] * Det3_345_045 +
1595 m[A14] * Det3_345_015 - m[A15] * Det3_345_014;
1596 G4double Det4_1345_0234 = m[A10] * Det3_345_234 - m[A12] * Det3_345_034 +
1597 m[A13] * Det3_345_024 - m[A14] * Det3_345_023;
1598 G4double Det4_1345_0235 = m[A10] * Det3_345_235 - m[A12] * Det3_345_035 +
1599 m[A13] * Det3_345_025 - m[A15] * Det3_345_023;
1600 G4double Det4_1345_0245 = m[A10] * Det3_345_245 - m[A12] * Det3_345_045 +
1601 m[A14] * Det3_345_025 - m[A15] * Det3_345_024;
1602 G4double Det4_1345_0345 = m[A10] * Det3_345_345 - m[A13] * Det3_345_045 +
1603 m[A14] * Det3_345_035 - m[A15] * Det3_345_034;
1604 G4double Det4_1345_1234 = m[A11] * Det3_345_234 - m[A12] * Det3_345_134 +
1605 m[A13] * Det3_345_124 - m[A14] * Det3_345_123;
1606 G4double Det4_1345_1235 = m[A11] * Det3_345_235 - m[A12] * Det3_345_135 +
1607 m[A13] * Det3_345_125 - m[A15] * Det3_345_123;
1608 G4double Det4_1345_1245 = m[A11] * Det3_345_245 - m[A12] * Det3_345_145 +
1609 m[A14] * Det3_345_125 - m[A15] * Det3_345_124;
1610 G4double Det4_1345_1345 = m[A11] * Det3_345_345 - m[A13] * Det3_345_145 +
1611 m[A14] * Det3_345_135 - m[A15] * Det3_345_134;
1612 G4double Det4_1345_2345 = m[A12] * Det3_345_345 - m[A13] * Det3_345_245 +
1613 m[A14] * Det3_345_235 - m[A15] * Det3_345_234;
1614 G4double Det4_2345_0123 = m[A20] * Det3_345_123 - m[A21] * Det3_345_023 +
1615 m[A22] * Det3_345_013 - m[A23] * Det3_345_012;
1616 G4double Det4_2345_0124 = m[A20] * Det3_345_124 - m[A21] * Det3_345_024 +
1617 m[A22] * Det3_345_014 - m[A24] * Det3_345_012;
1618 G4double Det4_2345_0125 = m[A20] * Det3_345_125 - m[A21] * Det3_345_025 +
1619 m[A22] * Det3_345_015 - m[A25] * Det3_345_012;
1620 G4double Det4_2345_0134 = m[A20] * Det3_345_134 - m[A21] * Det3_345_034 +
1621 m[A23] * Det3_345_014 - m[A24] * Det3_345_013;
1622 G4double Det4_2345_0135 = m[A20] * Det3_345_135 - m[A21] * Det3_345_035 +
1623 m[A23] * Det3_345_015 - m[A25] * Det3_345_013;
1624 G4double Det4_2345_0145 = m[A20] * Det3_345_145 - m[A21] * Det3_345_045 +
1625 m[A24] * Det3_345_015 - m[A25] * Det3_345_014;
1626 G4double Det4_2345_0234 = m[A20] * Det3_345_234 - m[A22] * Det3_345_034 +
1627 m[A23] * Det3_345_024 - m[A24] * Det3_345_023;
1628 G4double Det4_2345_0235 = m[A20] * Det3_345_235 - m[A22] * Det3_345_035 +
1629 m[A23] * Det3_345_025 - m[A25] * Det3_345_023;
1630 G4double Det4_2345_0245 = m[A20] * Det3_345_245 - m[A22] * Det3_345_045 +
1631 m[A24] * Det3_345_025 - m[A25] * Det3_345_024;
1632 G4double Det4_2345_0345 = m[A20] * Det3_345_345 - m[A23] * Det3_345_045 +
1633 m[A24] * Det3_345_035 - m[A25] * Det3_345_034;
1634 G4double Det4_2345_1234 = m[A21] * Det3_345_234 - m[A22] * Det3_345_134 +
1635 m[A23] * Det3_345_124 - m[A24] * Det3_345_123;
1636 G4double Det4_2345_1235 = m[A21] * Det3_345_235 - m[A22] * Det3_345_135 +
1637 m[A23] * Det3_345_125 - m[A25] * Det3_345_123;
1638 G4double Det4_2345_1245 = m[A21] * Det3_345_245 - m[A22] * Det3_345_145 +
1639 m[A24] * Det3_345_125 - m[A25] * Det3_345_124;
1640 G4double Det4_2345_1345 = m[A21] * Det3_345_345 - m[A23] * Det3_345_145 +
1641 m[A24] * Det3_345_135 - m[A25] * Det3_345_134;
1642 G4double Det4_2345_2345 = m[A22] * Det3_345_345 - m[A23] * Det3_345_245 +
1643 m[A24] * Det3_345_235 - m[A25] * Det3_345_234;
1644
1645 // Find all NECESSARY 5x5 dets: (36 of them)
1646
1647 G4double Det5_01234_01234 =
1648 m[A00] * Det4_1234_1234 - m[A01] * Det4_1234_0234 +
1649 m[A02] * Det4_1234_0134 - m[A03] * Det4_1234_0124 + m[A04] * Det4_1234_0123;
1650 G4double Det5_01234_01235 =
1651 m[A00] * Det4_1234_1235 - m[A01] * Det4_1234_0235 +
1652 m[A02] * Det4_1234_0135 - m[A03] * Det4_1234_0125 + m[A05] * Det4_1234_0123;
1653 G4double Det5_01234_01245 =
1654 m[A00] * Det4_1234_1245 - m[A01] * Det4_1234_0245 +
1655 m[A02] * Det4_1234_0145 - m[A04] * Det4_1234_0125 + m[A05] * Det4_1234_0124;
1656 G4double Det5_01234_01345 =
1657 m[A00] * Det4_1234_1345 - m[A01] * Det4_1234_0345 +
1658 m[A03] * Det4_1234_0145 - m[A04] * Det4_1234_0135 + m[A05] * Det4_1234_0134;
1659 G4double Det5_01234_02345 =
1660 m[A00] * Det4_1234_2345 - m[A02] * Det4_1234_0345 +
1661 m[A03] * Det4_1234_0245 - m[A04] * Det4_1234_0235 + m[A05] * Det4_1234_0234;
1662 G4double Det5_01234_12345 =
1663 m[A01] * Det4_1234_2345 - m[A02] * Det4_1234_1345 +
1664 m[A03] * Det4_1234_1245 - m[A04] * Det4_1234_1235 + m[A05] * Det4_1234_1234;
1665 G4double Det5_01235_01234 =
1666 m[A00] * Det4_1235_1234 - m[A01] * Det4_1235_0234 +
1667 m[A02] * Det4_1235_0134 - m[A03] * Det4_1235_0124 + m[A04] * Det4_1235_0123;
1668 G4double Det5_01235_01235 =
1669 m[A00] * Det4_1235_1235 - m[A01] * Det4_1235_0235 +
1670 m[A02] * Det4_1235_0135 - m[A03] * Det4_1235_0125 + m[A05] * Det4_1235_0123;
1671 G4double Det5_01235_01245 =
1672 m[A00] * Det4_1235_1245 - m[A01] * Det4_1235_0245 +
1673 m[A02] * Det4_1235_0145 - m[A04] * Det4_1235_0125 + m[A05] * Det4_1235_0124;
1674 G4double Det5_01235_01345 =
1675 m[A00] * Det4_1235_1345 - m[A01] * Det4_1235_0345 +
1676 m[A03] * Det4_1235_0145 - m[A04] * Det4_1235_0135 + m[A05] * Det4_1235_0134;
1677 G4double Det5_01235_02345 =
1678 m[A00] * Det4_1235_2345 - m[A02] * Det4_1235_0345 +
1679 m[A03] * Det4_1235_0245 - m[A04] * Det4_1235_0235 + m[A05] * Det4_1235_0234;
1680 G4double Det5_01235_12345 =
1681 m[A01] * Det4_1235_2345 - m[A02] * Det4_1235_1345 +
1682 m[A03] * Det4_1235_1245 - m[A04] * Det4_1235_1235 + m[A05] * Det4_1235_1234;
1683 G4double Det5_01245_01234 =
1684 m[A00] * Det4_1245_1234 - m[A01] * Det4_1245_0234 +
1685 m[A02] * Det4_1245_0134 - m[A03] * Det4_1245_0124 + m[A04] * Det4_1245_0123;
1686 G4double Det5_01245_01235 =
1687 m[A00] * Det4_1245_1235 - m[A01] * Det4_1245_0235 +
1688 m[A02] * Det4_1245_0135 - m[A03] * Det4_1245_0125 + m[A05] * Det4_1245_0123;
1689 G4double Det5_01245_01245 =
1690 m[A00] * Det4_1245_1245 - m[A01] * Det4_1245_0245 +
1691 m[A02] * Det4_1245_0145 - m[A04] * Det4_1245_0125 + m[A05] * Det4_1245_0124;
1692 G4double Det5_01245_01345 =
1693 m[A00] * Det4_1245_1345 - m[A01] * Det4_1245_0345 +
1694 m[A03] * Det4_1245_0145 - m[A04] * Det4_1245_0135 + m[A05] * Det4_1245_0134;
1695 G4double Det5_01245_02345 =
1696 m[A00] * Det4_1245_2345 - m[A02] * Det4_1245_0345 +
1697 m[A03] * Det4_1245_0245 - m[A04] * Det4_1245_0235 + m[A05] * Det4_1245_0234;
1698 G4double Det5_01245_12345 =
1699 m[A01] * Det4_1245_2345 - m[A02] * Det4_1245_1345 +
1700 m[A03] * Det4_1245_1245 - m[A04] * Det4_1245_1235 + m[A05] * Det4_1245_1234;
1701 G4double Det5_01345_01234 =
1702 m[A00] * Det4_1345_1234 - m[A01] * Det4_1345_0234 +
1703 m[A02] * Det4_1345_0134 - m[A03] * Det4_1345_0124 + m[A04] * Det4_1345_0123;
1704 G4double Det5_01345_01235 =
1705 m[A00] * Det4_1345_1235 - m[A01] * Det4_1345_0235 +
1706 m[A02] * Det4_1345_0135 - m[A03] * Det4_1345_0125 + m[A05] * Det4_1345_0123;
1707 G4double Det5_01345_01245 =
1708 m[A00] * Det4_1345_1245 - m[A01] * Det4_1345_0245 +
1709 m[A02] * Det4_1345_0145 - m[A04] * Det4_1345_0125 + m[A05] * Det4_1345_0124;
1710 G4double Det5_01345_01345 =
1711 m[A00] * Det4_1345_1345 - m[A01] * Det4_1345_0345 +
1712 m[A03] * Det4_1345_0145 - m[A04] * Det4_1345_0135 + m[A05] * Det4_1345_0134;
1713 G4double Det5_01345_02345 =
1714 m[A00] * Det4_1345_2345 - m[A02] * Det4_1345_0345 +
1715 m[A03] * Det4_1345_0245 - m[A04] * Det4_1345_0235 + m[A05] * Det4_1345_0234;
1716 G4double Det5_01345_12345 =
1717 m[A01] * Det4_1345_2345 - m[A02] * Det4_1345_1345 +
1718 m[A03] * Det4_1345_1245 - m[A04] * Det4_1345_1235 +
1719 m[A05] * Det4_1345_1234; //
1720 G4double Det5_02345_01234 =
1721 m[A00] * Det4_2345_1234 - m[A01] * Det4_2345_0234 +
1722 m[A02] * Det4_2345_0134 - m[A03] * Det4_2345_0124 + m[A04] * Det4_2345_0123;
1723 G4double Det5_02345_01235 =
1724 m[A00] * Det4_2345_1235 - m[A01] * Det4_2345_0235 +
1725 m[A02] * Det4_2345_0135 - m[A03] * Det4_2345_0125 + m[A05] * Det4_2345_0123;
1726 G4double Det5_02345_01245 =
1727 m[A00] * Det4_2345_1245 - m[A01] * Det4_2345_0245 +
1728 m[A02] * Det4_2345_0145 - m[A04] * Det4_2345_0125 + m[A05] * Det4_2345_0124;
1729 G4double Det5_02345_01345 =
1730 m[A00] * Det4_2345_1345 - m[A01] * Det4_2345_0345 +
1731 m[A03] * Det4_2345_0145 - m[A04] * Det4_2345_0135 + m[A05] * Det4_2345_0134;
1732 G4double Det5_02345_02345 =
1733 m[A00] * Det4_2345_2345 - m[A02] * Det4_2345_0345 +
1734 m[A03] * Det4_2345_0245 - m[A04] * Det4_2345_0235 + m[A05] * Det4_2345_0234;
1735 G4double Det5_02345_12345 =
1736 m[A01] * Det4_2345_2345 - m[A02] * Det4_2345_1345 +
1737 m[A03] * Det4_2345_1245 - m[A04] * Det4_2345_1235 + m[A05] * Det4_2345_1234;
1738 G4double Det5_12345_01234 =
1739 m[A10] * Det4_2345_1234 - m[A11] * Det4_2345_0234 +
1740 m[A12] * Det4_2345_0134 - m[A13] * Det4_2345_0124 + m[A14] * Det4_2345_0123;
1741 G4double Det5_12345_01235 =
1742 m[A10] * Det4_2345_1235 - m[A11] * Det4_2345_0235 +
1743 m[A12] * Det4_2345_0135 - m[A13] * Det4_2345_0125 + m[A15] * Det4_2345_0123;
1744 G4double Det5_12345_01245 =
1745 m[A10] * Det4_2345_1245 - m[A11] * Det4_2345_0245 +
1746 m[A12] * Det4_2345_0145 - m[A14] * Det4_2345_0125 + m[A15] * Det4_2345_0124;
1747 G4double Det5_12345_01345 =
1748 m[A10] * Det4_2345_1345 - m[A11] * Det4_2345_0345 +
1749 m[A13] * Det4_2345_0145 - m[A14] * Det4_2345_0135 + m[A15] * Det4_2345_0134;
1750 G4double Det5_12345_02345 =
1751 m[A10] * Det4_2345_2345 - m[A12] * Det4_2345_0345 +
1752 m[A13] * Det4_2345_0245 - m[A14] * Det4_2345_0235 + m[A15] * Det4_2345_0234;
1753 G4double Det5_12345_12345 =
1754 m[A11] * Det4_2345_2345 - m[A12] * Det4_2345_1345 +
1755 m[A13] * Det4_2345_1245 - m[A14] * Det4_2345_1235 + m[A15] * Det4_2345_1234;
1756
1757 // Find the determinant
1758
1759 G4double det = m[A00] * Det5_12345_12345 - m[A01] * Det5_12345_02345 +
1760 m[A02] * Det5_12345_01345 - m[A03] * Det5_12345_01245 +
1761 m[A04] * Det5_12345_01235 - m[A05] * Det5_12345_01234;
1762
1763 if(det == 0)
1764 {
1765 ifail = 1;
1766 return;
1767 }
1768
1769 G4double oneOverDet = 1.0 / det;
1770 G4double mn1OverDet = -oneOverDet;
1771
1772 m[A00] = Det5_12345_12345 * oneOverDet;
1773 m[A01] = Det5_02345_12345 * mn1OverDet;
1774 m[A02] = Det5_01345_12345 * oneOverDet;
1775 m[A03] = Det5_01245_12345 * mn1OverDet;
1776 m[A04] = Det5_01235_12345 * oneOverDet;
1777 m[A05] = Det5_01234_12345 * mn1OverDet;
1778
1779 m[A10] = Det5_12345_02345 * mn1OverDet;
1780 m[A11] = Det5_02345_02345 * oneOverDet;
1781 m[A12] = Det5_01345_02345 * mn1OverDet;
1782 m[A13] = Det5_01245_02345 * oneOverDet;
1783 m[A14] = Det5_01235_02345 * mn1OverDet;
1784 m[A15] = Det5_01234_02345 * oneOverDet;
1785
1786 m[A20] = Det5_12345_01345 * oneOverDet;
1787 m[A21] = Det5_02345_01345 * mn1OverDet;
1788 m[A22] = Det5_01345_01345 * oneOverDet;
1789 m[A23] = Det5_01245_01345 * mn1OverDet;
1790 m[A24] = Det5_01235_01345 * oneOverDet;
1791 m[A25] = Det5_01234_01345 * mn1OverDet;
1792
1793 m[A30] = Det5_12345_01245 * mn1OverDet;
1794 m[A31] = Det5_02345_01245 * oneOverDet;
1795 m[A32] = Det5_01345_01245 * mn1OverDet;
1796 m[A33] = Det5_01245_01245 * oneOverDet;
1797 m[A34] = Det5_01235_01245 * mn1OverDet;
1798 m[A35] = Det5_01234_01245 * oneOverDet;
1799
1800 m[A40] = Det5_12345_01235 * oneOverDet;
1801 m[A41] = Det5_02345_01235 * mn1OverDet;
1802 m[A42] = Det5_01345_01235 * oneOverDet;
1803 m[A43] = Det5_01245_01235 * mn1OverDet;
1804 m[A44] = Det5_01235_01235 * oneOverDet;
1805 m[A45] = Det5_01234_01235 * mn1OverDet;
1806
1807 m[A50] = Det5_12345_01234 * mn1OverDet;
1808 m[A51] = Det5_02345_01234 * oneOverDet;
1809 m[A52] = Det5_01345_01234 * mn1OverDet;
1810 m[A53] = Det5_01245_01234 * oneOverDet;
1811 m[A54] = Det5_01235_01234 * mn1OverDet;
1812 m[A55] = Det5_01234_01234 * oneOverDet;
1813
1814 return;
1815}
A41
#define A41
Definition: G4ErrorMatrix.cc:896
A20
#define A20
Definition: G4ErrorMatrix.cc:881
A01
#define A01
Definition: G4ErrorMatrix.cc:868
A53
#define A53
Definition: G4ErrorMatrix.cc:905
A23
#define A23
Definition: G4ErrorMatrix.cc:884
A45
#define A45
Definition: G4ErrorMatrix.cc:900
A13
#define A13
Definition: G4ErrorMatrix.cc:877
A00
#define A00
Definition: G4ErrorMatrix.cc:867
A54
#define A54
Definition: G4ErrorMatrix.cc:906
A40
#define A40
Definition: G4ErrorMatrix.cc:895
A25
#define A25
Definition: G4ErrorMatrix.cc:886
A02
#define A02
Definition: G4ErrorMatrix.cc:869
A24
#define A24
Definition: G4ErrorMatrix.cc:885
A22
#define A22
Definition: G4ErrorMatrix.cc:883
A04
#define A04
Definition: G4ErrorMatrix.cc:871
A30
#define A30
Definition: G4ErrorMatrix.cc:888
A12
#define A12
Definition: G4ErrorMatrix.cc:876
A55
#define A55
Definition: G4ErrorMatrix.cc:907
A35
#define A35
Definition: G4ErrorMatrix.cc:893
A50
#define A50
Definition: G4ErrorMatrix.cc:902
A03
#define A03
Definition: G4ErrorMatrix.cc:870
A14
#define A14
Definition: G4ErrorMatrix.cc:878
A51
#define A51
Definition: G4ErrorMatrix.cc:903
A21
#define A21
Definition: G4ErrorMatrix.cc:882
A11
#define A11
Definition: G4ErrorMatrix.cc:875
A10
#define A10
Definition: G4ErrorMatrix.cc:874
A44
#define A44
Definition: G4ErrorMatrix.cc:899
A05
#define A05
Definition: G4ErrorMatrix.cc:872
A32
#define A32
Definition: G4ErrorMatrix.cc:890
A31
#define A31
Definition: G4ErrorMatrix.cc:889
A33
#define A33
Definition: G4ErrorMatrix.cc:891
A42
#define A42
Definition: G4ErrorMatrix.cc:897
A52
#define A52
Definition: G4ErrorMatrix.cc:904
A15
#define A15
Definition: G4ErrorMatrix.cc:879
A34
#define A34
Definition: G4ErrorMatrix.cc:892
A43
#define A43
Definition: G4ErrorMatrix.cc:898

References A00, A01, A02, A03, A04, A05, A10, A11, A12, A13, A14, A15, A20, A21, A22, A23, A24, A25, A30, A31, A32, A33, A34, A35, A40, A41, A42, A43, A44, A45, A50, A51, A52, A53, A54, A55, and m.

Referenced by invert().

◆ num_col()

virtual G4int G4ErrorMatrix::num_col ( ) const
inlinevirtual

Referenced by apply(), dfinv_matrix(), dsum(), operator+(), operator+=(), operator-(), operator-=(), operator<<(), G4ErrorSymMatrix::similarity(), G4ErrorSymMatrix::similarityT(), and sub().

◆ num_row()

virtual G4int G4ErrorMatrix::num_row ( ) const
inlinevirtual

Referenced by apply(), dfinv_matrix(), dsum(), operator+(), operator+=(), operator-(), operator-=(), operator<<(), operator=(), G4ErrorSymMatrix::similarity(), G4ErrorSymMatrix::similarityT(), and sub().

◆ num_size()

virtual G4int G4ErrorMatrix::num_size ( ) const
inlineprotectedvirtual

◆ operator()() [1/2]

virtual G4double & G4ErrorMatrix::operator() ( G4int  row,
G4int  col 
)
inlinevirtual

◆ operator()() [2/2]

virtual const G4double & G4ErrorMatrix::operator() ( G4int  row,
G4int  col 
) const
inlinevirtual

◆ operator*=()

G4ErrorMatrix & G4ErrorMatrix::operator*= ( G4double  t)

Definition at line 345 of file G4ErrorMatrix.cc.

346{
347 SIMPLE_UOP(*=)
348 return (*this);
349}
SIMPLE_UOP
#define SIMPLE_UOP(OPER)
Definition: G4ErrorMatrix.cc:41

References SIMPLE_UOP.

◆ operator+=() [1/2]

G4ErrorMatrix & G4ErrorMatrix::operator+= ( const G4ErrorMatrix &  m2)

Definition at line 325 of file G4ErrorMatrix.cc.

326{
327 CHK_DIM_2(num_row(), mat2.num_row(), num_col(), mat2.num_col(), +=);
328 SIMPLE_BOP(+=)
329 return (*this);
330}
CHK_DIM_2
#define CHK_DIM_2(r1, r2, c1, c2, fun)
Definition: G4ErrorMatrix.cc:64
SIMPLE_BOP
#define SIMPLE_BOP(OPER)
Definition: G4ErrorMatrix.cc:47

References CHK_DIM_2, num_col(), num_row(), and SIMPLE_BOP.

◆ operator+=() [2/2]

G4ErrorMatrix & G4ErrorMatrix::operator+= ( const G4ErrorSymMatrix &  m2)

Definition at line 436 of file G4ErrorSymMatrix.cc.

437{
438 CHK_DIM_2(num_row(), mat2.num_row(), num_col(), mat2.num_col(), +=);
439 G4int n = num_col();
440 G4ErrorMatrixConstIter sjk = mat2.m.begin();
441 G4ErrorMatrixIter m1j = m.begin();
442 G4ErrorMatrixIter mj = m.begin();
443 // j >= k
444 for(G4int j = 1; j <= num_row(); j++)
445 {
446 G4ErrorMatrixIter mjk = mj;
447 G4ErrorMatrixIter mkj = m1j;
448 for(G4int k = 1; k <= j; k++)
449 {
450 *(mjk++) += *sjk;
451 if(j != k)
452 *mkj += *sjk;
453 sjk++;
454 mkj += n;
455 }
456 mj += n;
457 m1j++;
458 }
459 return (*this);
460}
CHK_DIM_2
#define CHK_DIM_2(r1, r2, c1, c2, fun)
Definition: G4ErrorSymMatrix.cc:61

References CHK_DIM_2, m, G4ErrorSymMatrix::m, CLHEP::detail::n, num_col(), G4ErrorSymMatrix::num_col(), num_row(), and G4ErrorSymMatrix::num_row().

◆ operator-()

G4ErrorMatrix G4ErrorMatrix::operator- ( ) const

Definition at line 234 of file G4ErrorMatrix.cc.

235{
236 G4ErrorMatrix mat2(nrow, ncol);
237 G4ErrorMatrixConstIter a = m.begin();
238 G4ErrorMatrixIter b = mat2.m.begin();
239 G4ErrorMatrixConstIter e = m.end();
240 for(; a < e; a++, b++)
241 (*b) = -(*a);
242 return mat2;
243}

References m, ncol, and nrow.

◆ operator-=() [1/2]

G4ErrorMatrix & G4ErrorMatrix::operator-= ( const G4ErrorMatrix &  m2)

Definition at line 332 of file G4ErrorMatrix.cc.

333{
334 CHK_DIM_2(num_row(), mat2.num_row(), num_col(), mat2.num_col(), -=);
335 SIMPLE_BOP(-=)
336 return (*this);
337}

References CHK_DIM_2, num_col(), num_row(), and SIMPLE_BOP.

◆ operator-=() [2/2]

G4ErrorMatrix & G4ErrorMatrix::operator-= ( const G4ErrorSymMatrix &  m2)

Definition at line 469 of file G4ErrorSymMatrix.cc.

470{
471 CHK_DIM_2(num_row(), mat2.num_row(), num_col(), mat2.num_col(), -=);
472 G4int n = num_col();
473 G4ErrorMatrixConstIter sjk = mat2.m.begin();
474 G4ErrorMatrixIter m1j = m.begin();
475 G4ErrorMatrixIter mj = m.begin();
476 // j >= k
477 for(G4int j = 1; j <= num_row(); j++)
478 {
479 G4ErrorMatrixIter mjk = mj;
480 G4ErrorMatrixIter mkj = m1j;
481 for(G4int k = 1; k <= j; k++)
482 {
483 *(mjk++) -= *sjk;
484 if(j != k)
485 *mkj -= *sjk;
486 sjk++;
487 mkj += n;
488 }
489 mj += n;
490 m1j++;
491 }
492 return (*this);
493}

References CHK_DIM_2, m, G4ErrorSymMatrix::m, CLHEP::detail::n, num_col(), G4ErrorSymMatrix::num_col(), num_row(), and G4ErrorSymMatrix::num_row().

◆ operator/=()

G4ErrorMatrix & G4ErrorMatrix::operator/= ( G4double  t)

Definition at line 339 of file G4ErrorMatrix.cc.

340{
341 SIMPLE_UOP(/=)
342 return (*this);
343}

References SIMPLE_UOP.

◆ operator=() [1/3]

G4ErrorMatrix & G4ErrorMatrix::operator= ( const G4ErrorMatrix &  m2)

Definition at line 351 of file G4ErrorMatrix.cc.

352{
353 if(&mat1 == this)
354 {
355 return *this;
356 }
357
358 if(mat1.nrow * mat1.ncol != size) //??fixme?? mat1.size != size
359 {
360 size = mat1.nrow * mat1.ncol;
361 m.resize(size); //??fixme?? if (size < mat1.size) m.resize(mat1.size);
362 }
363 nrow = mat1.nrow;
364 ncol = mat1.ncol;
365 m = mat1.m;
366 return (*this);
367}

References m, ncol, nrow, and size.

◆ operator=() [2/3]

G4ErrorMatrix & G4ErrorMatrix::operator= ( const G4ErrorSymMatrix &  m2)

Definition at line 514 of file G4ErrorSymMatrix.cc.

515{
516 if(mat1.nrow * mat1.nrow != size)
517 {
518 size = mat1.nrow * mat1.nrow;
519 m.resize(size);
520 }
521 nrow = mat1.nrow;
522 ncol = mat1.nrow;
523 G4int n = ncol;
524 G4ErrorMatrixConstIter sjk = mat1.m.begin();
525 G4ErrorMatrixIter m1j = m.begin();
526 G4ErrorMatrixIter mj = m.begin();
527 // j >= k
528 for(G4int j = 1; j <= num_row(); j++)
529 {
530 G4ErrorMatrixIter mjk = mj;
531 G4ErrorMatrixIter mkj = m1j;
532 for(G4int k = 1; k <= j; k++)
533 {
534 *(mjk++) = *sjk;
535 if(j != k)
536 *mkj = *sjk;
537 sjk++;
538 mkj += n;
539 }
540 mj += n;
541 m1j++;
542 }
543 return (*this);
544}

References m, G4ErrorSymMatrix::m, CLHEP::detail::n, ncol, nrow, G4ErrorSymMatrix::nrow, num_row(), and size.

◆ operator=() [3/3]

G4ErrorMatrix & G4ErrorMatrix::operator= ( G4ErrorMatrix &&  )
default

◆ operator[]() [1/2]

G4ErrorMatrix_row G4ErrorMatrix::operator[] ( G4int  )
inline

◆ operator[]() [2/2]

const G4ErrorMatrix_row_const G4ErrorMatrix::operator[] ( G4int  ) const
inline

◆ sub() [1/2]

G4ErrorMatrix G4ErrorMatrix::sub ( G4int  min_row,
G4int  max_row,
G4int  min_col,
G4int  max_col 
) const

Definition at line 170 of file G4ErrorMatrix.cc.

172{
173 G4ErrorMatrix mret(max_row - min_row + 1, max_col - min_col + 1);
174 if(max_row > num_row() || max_col > num_col())
175 {
176 error("G4ErrorMatrix::sub: Index out of range");
177 }
178 G4ErrorMatrixIter a = mret.m.begin();
179 G4int nc = num_col();
180 G4ErrorMatrixConstIter b1 = m.begin() + (min_row - 1) * nc + min_col - 1;
181
182 for(G4int irow = 1; irow <= mret.num_row(); irow++)
183 {
184 G4ErrorMatrixConstIter brc = b1;
185 for(G4int icol = 1; icol <= mret.num_col(); icol++)
186 {
187 *(a++) = *(brc++);
188 }
189 b1 += nc;
190 }
191 return mret;
192}

References error(), m, num_col(), and num_row().

Referenced by dsum().

◆ sub() [2/2]

void G4ErrorMatrix::sub ( G4int  row,
G4int  col,
const G4ErrorMatrix &  m1 
)

Definition at line 194 of file G4ErrorMatrix.cc.

195{
196 if(row < 1 || row + mat1.num_row() - 1 > num_row() || col < 1 ||
197 col + mat1.num_col() - 1 > num_col())
198 {
199 error("G4ErrorMatrix::sub: Index out of range");
200 }
201 G4ErrorMatrixConstIter a = mat1.m.begin();
202 G4int nc = num_col();
203 G4ErrorMatrixIter b1 = m.begin() + (row - 1) * nc + col - 1;
204
205 for(G4int irow = 1; irow <= mat1.num_row(); irow++)
206 {
207 G4ErrorMatrixIter brc = b1;
208 for(G4int icol = 1; icol <= mat1.num_col(); icol++)
209 {
210 *(brc++) = *(a++);
211 }
212 b1 += nc;
213 }
214}

References error(), m, num_col(), and num_row().

◆ T()

G4ErrorMatrix G4ErrorMatrix::T ( ) const

Definition at line 399 of file G4ErrorMatrix.cc.

400{
401 G4ErrorMatrix mret(ncol, nrow);
402 G4ErrorMatrixConstIter pl = m.end();
403 G4ErrorMatrixConstIter pme = m.begin();
404 G4ErrorMatrixIter pt = mret.m.begin();
405 G4ErrorMatrixIter ptl = mret.m.end();
406 for(; pme < pl; pme++, pt += nrow)
407 {
408 if(pt >= ptl)
409 {
410 pt -= (size - 1);
411 }
412 (*pt) = (*pme);
413 }
414 return mret;
415}

References m, ncol, nrow, and size.

Referenced by G4ErrorFreeTrajState::PropagateError().

◆ trace()

G4double G4ErrorMatrix::trace ( ) const

Definition at line 845 of file G4ErrorMatrix.cc.

846{
847 G4double t = 0.0;
848 for(G4ErrorMatrixConstIter d = m.begin(); d < m.end(); d += (ncol + 1))
849 {
850 t += *d;
851 }
852 return t;
853}

References m, and ncol.

Friends And Related Function Documentation

◆ back_solve

void back_solve ( const G4ErrorMatrix &  R,
G4ErrorMatrix *  b 
)
friend

◆ col_givens

void col_givens ( G4ErrorMatrix *  A,
G4double  c,
G4double  s,
G4int  k1,
G4int  k2,
G4int  rowmin,
G4int  rowmax 
)
friend

◆ col_house

void col_house ( G4ErrorMatrix *  ,
const G4ErrorMatrix &  ,
G4double  ,
G4int  ,
G4int  ,
G4int  ,
G4int   
)
friend

◆ G4ErrorMatrix_row

friend class G4ErrorMatrix_row
friend

Definition at line 174 of file G4ErrorMatrix.hh.

◆ G4ErrorMatrix_row_const

friend class G4ErrorMatrix_row_const
friend

Definition at line 175 of file G4ErrorMatrix.hh.

◆ G4ErrorSymMatrix

friend class G4ErrorSymMatrix
friend

Definition at line 176 of file G4ErrorMatrix.hh.

◆ house_with_update [1/2]

void house_with_update ( G4ErrorMatrix *  a,
G4ErrorMatrix *  v,
G4int  row,
G4int  col 
)
friend

◆ house_with_update [2/2]

void house_with_update ( G4ErrorMatrix *  a,
G4int  row,
G4int  col 
)
friend

◆ house_with_update2

void house_with_update2 ( G4ErrorSymMatrix *  a,
G4ErrorMatrix *  v,
G4int  row,
G4int  col 
)
friend

◆ operator* [1/4]

G4ErrorMatrix operator* ( const G4ErrorMatrix &  m1,
const G4ErrorMatrix &  m2 
)
friend

Definition at line 291 of file G4ErrorMatrix.cc.

292{
293 // initialize matrix to 0.0
294 G4ErrorMatrix mret(mat1.nrow, mat2.ncol, 0);
295 CHK_DIM_1(mat1.ncol, mat2.nrow, *);
296
297 G4int m1cols = mat1.ncol;
298 G4int m2cols = mat2.ncol;
299
300 for(G4int i = 0; i < mat1.nrow; i++)
301 {
302 for(G4int j = 0; j < m1cols; j++)
303 {
304 G4double temp = mat1.m[i * m1cols + j];
305 G4ErrorMatrixIter pt = mret.m.begin() + i * m2cols;
306
307 // Loop over k (the column index in matrix mat2)
308 G4ErrorMatrixConstIter pb = mat2.m.begin() + m2cols * j;
309 const G4ErrorMatrixConstIter pblast = pb + m2cols;
310 while(pb < pblast) // Loop checking, 06.08.2015, G.Cosmo
311 {
312 (*pt) += temp * (*pb);
313 pb++;
314 pt++;
315 }
316 }
317 }
318 return mret;
319}
CHK_DIM_1
#define CHK_DIM_1(c1, r2, fun)
Definition: G4ErrorMatrix.cc:70

◆ operator* [2/4]

G4ErrorMatrix operator* ( const G4ErrorMatrix &  m1,
const G4ErrorSymMatrix &  m2 
)
friend

Definition at line 302 of file G4ErrorSymMatrix.cc.

303{
304 G4ErrorMatrix mret(mat1.num_row(), mat2.num_col());
305 CHK_DIM_1(mat1.num_col(), mat2.num_row(), *);
306 G4ErrorMatrixConstIter mit1, mit2, sp, snp; // mit2=0
307 G4double temp;
308 G4ErrorMatrixIter mir = mret.m.begin();
309 for(mit1 = mat1.m.begin();
310 mit1 < mat1.m.begin() + mat1.num_row() * mat1.num_col(); mit1 = mit2)
311 {
312 snp = mat2.m.begin();
313 for(int step = 1; step <= mat2.num_row(); ++step)
314 {
315 mit2 = mit1;
316 sp = snp;
317 snp += step;
318 temp = 0;
319 while(sp < snp) // Loop checking, 06.08.2015, G.Cosmo
320 {
321 temp += *(sp++) * (*(mit2++));
322 }
323 if(step < mat2.num_row())
324 { // only if we aren't on the last row
325 sp += step - 1;
326 for(int stept = step + 1; stept <= mat2.num_row(); stept++)
327 {
328 temp += *sp * (*(mit2++));
329 if(stept < mat2.num_row())
330 sp += stept;
331 }
332 } // if(step
333 *(mir++) = temp;
334 } // for(step
335 } // for(mit1
336 return mret;
337}
CHK_DIM_1
#define CHK_DIM_1(c1, r2, fun)
Definition: G4ErrorSymMatrix.cc:67

◆ operator* [3/4]

G4ErrorMatrix operator* ( const G4ErrorSymMatrix &  m1,
const G4ErrorMatrix &  m2 
)
friend

Definition at line 339 of file G4ErrorSymMatrix.cc.

340{
341 G4ErrorMatrix mret(mat1.num_row(), mat2.num_col());
342 CHK_DIM_1(mat1.num_col(), mat2.num_row(), *);
343 G4int step, stept;
344 G4ErrorMatrixConstIter mit1, mit2, sp, snp;
345 G4double temp;
346 G4ErrorMatrixIter mir = mret.m.begin();
347 for(step = 1, snp = mat1.m.begin(); step <= mat1.num_row(); snp += step++)
348 {
349 for(mit1 = mat2.m.begin(); mit1 < mat2.m.begin() + mat2.num_col(); mit1++)
350 {
351 mit2 = mit1;
352 sp = snp;
353 temp = 0;
354 while(sp < snp + step) // Loop checking, 06.08.2015, G.Cosmo
355 {
356 temp += *mit2 * (*(sp++));
357 mit2 += mat2.num_col();
358 }
359 sp += step - 1;
360 for(stept = step + 1; stept <= mat1.num_row(); stept++)
361 {
362 temp += *mit2 * (*sp);
363 mit2 += mat2.num_col();
364 sp += stept;
365 }
366 *(mir++) = temp;
367 }
368 }
369 return mret;
370}

◆ operator* [4/4]

G4ErrorMatrix operator* ( const G4ErrorSymMatrix &  m1,
const G4ErrorSymMatrix &  m2 
)
friend

Definition at line 372 of file G4ErrorSymMatrix.cc.

374{
375 G4ErrorMatrix mret(mat1.num_row(), mat1.num_row());
376 CHK_DIM_1(mat1.num_col(), mat2.num_row(), *);
377 G4int step1, stept1, step2, stept2;
378 G4ErrorMatrixConstIter snp1, sp1, snp2, sp2;
379 G4double temp;
380 G4ErrorMatrixIter mr = mret.m.begin();
381 for(step1 = 1, snp1 = mat1.m.begin(); step1 <= mat1.num_row();
382 snp1 += step1++)
383 {
384 for(step2 = 1, snp2 = mat2.m.begin(); step2 <= mat2.num_row();)
385 {
386 sp1 = snp1;
387 sp2 = snp2;
388 snp2 += step2;
389 temp = 0;
390 if(step1 < step2)
391 {
392 while(sp1 < snp1 + step1) // Loop checking, 06.08.2015, G.Cosmo
393 {
394 temp += (*(sp1++)) * (*(sp2++));
395 }
396 sp1 += step1 - 1;
397 for(stept1 = step1 + 1; stept1 != step2 + 1; sp1 += stept1++)
398 {
399 temp += (*sp1) * (*(sp2++));
400 }
401 sp2 += step2 - 1;
402 for(stept2 = ++step2; stept2 <= mat2.num_row();
403 sp1 += stept1++, sp2 += stept2++)
404 {
405 temp += (*sp1) * (*sp2);
406 }
407 }
408 else
409 {
410 while(sp2 < snp2) // Loop checking, 06.08.2015, G.Cosmo
411 {
412 temp += (*(sp1++)) * (*(sp2++));
413 }
414 sp2 += step2 - 1;
415 for(stept2 = ++step2; stept2 != step1 + 1; sp2 += stept2++)
416 {
417 temp += (*(sp1++)) * (*sp2);
418 }
419 sp1 += step1 - 1;
420 for(stept1 = step1 + 1; stept1 <= mat1.num_row();
421 sp1 += stept1++, sp2 += stept2++)
422 {
423 temp += (*sp1) * (*sp2);
424 }
425 }
426 *(mr++) = temp;
427 }
428 }
429 return mret;
430}

◆ operator+

G4ErrorMatrix operator+ ( const G4ErrorMatrix &  m1,
const G4ErrorMatrix &  m2 
)
friend

Definition at line 245 of file G4ErrorMatrix.cc.

246{
247 G4ErrorMatrix mret(mat1.nrow, mat1.ncol);
248 CHK_DIM_2(mat1.num_row(), mat2.num_row(), mat1.num_col(), mat2.num_col(), +);
249 SIMPLE_TOP(+)
250 return mret;
251}
SIMPLE_TOP
#define SIMPLE_TOP(OPER)
Definition: G4ErrorMatrix.cc:54

◆ operator-

G4ErrorMatrix operator- ( const G4ErrorMatrix &  m1,
const G4ErrorMatrix &  m2 
)
friend

Definition at line 257 of file G4ErrorMatrix.cc.

258{
259 G4ErrorMatrix mret(mat1.num_row(), mat1.num_col());
260 CHK_DIM_2(mat1.num_row(), mat2.num_row(), mat1.num_col(), mat2.num_col(), -);
261 SIMPLE_TOP(-)
262 return mret;
263}

◆ qr_solve

G4ErrorMatrix qr_solve ( G4ErrorMatrix *  ,
const G4ErrorMatrix &  b 
)
friend

◆ row_givens

void row_givens ( G4ErrorMatrix *  A,
G4double  c,
G4double  s,
G4int  k1,
G4int  k2,
G4int  colmin,
G4int  colmax 
)
friend

◆ row_house

void row_house ( G4ErrorMatrix *  ,
const G4ErrorMatrix &  ,
G4double  ,
G4int  ,
G4int  ,
G4int  ,
G4int   
)
friend

◆ tridiagonal

void tridiagonal ( G4ErrorSymMatrix *  a,
G4ErrorMatrix *  hsm 
)
friend

Field Documentation

◆ m

std::vector<G4double> G4ErrorMatrix::m
private

Definition at line 222 of file G4ErrorMatrix.hh.

Referenced by apply(), G4ErrorSymMatrix::assign(), dfact_matrix(), dfinv_matrix(), G4ErrorMatrix(), invert(), invertHaywood4(), invertHaywood5(), invertHaywood6(), operator+=(), operator-(), operator-=(), operator=(), G4ErrorSymMatrix::similarity(), G4ErrorSymMatrix::similarityT(), sub(), T(), and trace().

◆ ncol

G4int G4ErrorMatrix::ncol
private

Definition at line 224 of file G4ErrorMatrix.hh.

Referenced by determinant(), dfact_matrix(), G4ErrorMatrix(), invert(), operator-(), operator=(), T(), and trace().

◆ nrow

G4int G4ErrorMatrix::nrow
private

Definition at line 224 of file G4ErrorMatrix.hh.

Referenced by G4ErrorSymMatrix::assign(), determinant(), dfact_matrix(), G4ErrorMatrix(), invert(), operator-(), operator=(), and T().

◆ size

G4int G4ErrorMatrix::size
private

Definition at line 225 of file G4ErrorMatrix.hh.

Referenced by G4ErrorMatrix(), operator=(), and T().

The documentation for this class was generated from the following files:
Generated by doxygen 1.9.3