Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4ErrorSymMatrix.hh
이 파일의 문서화 페이지로 가기
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4ErrorSymMatrix.hh 108074 2017-12-19 15:35:08Z gcosmo $
27 //
28 // Class Description:
29 //
30 // Simplified version of CLHEP HepSymMatrix class.
31 
32 // History:
33 // - Imported from CLHEP and modified: P. Arce May 2007
34 // --------------------------------------------------------------------
35 
36 #ifndef G4ErrorSymMatrix_hh
37 #define G4ErrorSymMatrix_hh
38 
39 #include <vector>
40 #include "globals.hh"
41 
42 class G4ErrorMatrix;
43 
45 {
46  public: //with description
47 
48  inline G4ErrorSymMatrix();
49  // Default constructor. Gives 0x0 symmetric matrix.
50  // Another G4ErrorSymMatrix can be assigned to it.
51 
52  explicit G4ErrorSymMatrix(G4int p);
54  // Constructor. Gives p x p symmetric matrix.
55  // With a second argument, the matrix is initialized. 0 means a zero
56  // matrix, 1 means the identity matrix.
57 
59  // Copy constructor.
60 
61  // Constructor from DiagMatrix
62 
63  virtual ~G4ErrorSymMatrix();
64  // Destructor.
65 
66  inline G4int num_row() const;
67  inline G4int num_col() const;
68  // Returns number of rows/columns.
69 
70  const G4double & operator()(G4int row, G4int col) const;
71  G4double & operator()(G4int row, G4int col);
72  // Read and write a G4ErrorSymMatrix element.
73  // ** Note that indexing starts from (1,1). **
74 
75  const G4double & fast(G4int row, G4int col) const;
76  G4double & fast(G4int row, G4int col);
77  // fast element access.
78  // Must be row>=col;
79  // ** Note that indexing starts from (1,1). **
80 
81  void assign(const G4ErrorMatrix &m2);
82  // Assigns m2 to s, assuming m2 is a symmetric matrix.
83 
84  void assign(const G4ErrorSymMatrix &m2);
85  // Another form of assignment. For consistency.
86 
88  // Multiply a G4ErrorSymMatrix by a floating number.
89 
91  // Divide a G4ErrorSymMatrix by a floating number.
92 
95  // Add or subtract a G4ErrorSymMatrix.
96 
98  // Assignment operators. Notice that there is no G4ErrorSymMatrix = Matrix.
99 
101  // unary minus, ie. flip the sign of each element.
102 
103  G4ErrorSymMatrix T() const;
104  // Returns the transpose of a G4ErrorSymMatrix (which is itself).
105 
107  // Apply a function to all elements of the matrix.
108 
109  G4ErrorSymMatrix similarity(const G4ErrorMatrix &m1) const;
111  // Returns m1*s*m1.T().
112 
113  G4ErrorSymMatrix similarityT(const G4ErrorMatrix &m1) const;
114  // temporary. test of new similarity.
115  // Returns m1.T()*s*m1.
116 
117  G4ErrorSymMatrix sub(G4int min_row, G4int max_row) const;
118  // Returns a sub matrix of a G4ErrorSymMatrix.
119 
120  void sub(G4int row, const G4ErrorSymMatrix &m1);
121  // Sub matrix of this G4ErrorSymMatrix is replaced with m1.
122 
123  G4ErrorSymMatrix sub(G4int min_row, G4int max_row);
124  // SGI CC bug. I have to have both with/without const. I should not need
125  // one without const.
126 
127  inline G4ErrorSymMatrix inverse(G4int &ifail) const;
128  // Invert a Matrix. The matrix is not changed
129  // Returns 0 when successful, otherwise non-zero.
130 
131  void invert(G4int &ifail);
132  // Invert a Matrix.
133  // N.B. the contents of the matrix are replaced by the inverse.
134  // Returns ierr = 0 when successful, otherwise non-zero.
135  // This method has less overhead then inverse().
136 
137  G4double determinant() const;
138  // calculate the determinant of the matrix.
139 
140  G4double trace() const;
141  // calculate the trace of the matrix (sum of diagonal elements).
142 
144  {
145  public:
147  inline G4double & operator[](G4int);
148  private:
151  };
153  {
154  public:
156  inline const G4double & operator[](G4int) const;
157  private:
160  };
161  // helper class to implement m[i][j]
162 
165  // Read or write a matrix element.
166  // While it may not look like it, you simply do m[i][j] to get an
167  // element.
168  // ** Note that the indexing starts from [0][0]. **
169 
170  // Special-case inversions for 5x5 and 6x6 symmetric positive definite:
171  // These set ifail=0 and invert if the matrix was positive definite;
172  // otherwise ifail=1 and the matrix is left unaltered.
173 
174  void invertCholesky5 (G4int &ifail);
175  void invertCholesky6 (G4int &ifail);
176 
177  // Inversions for 5x5 and 6x6 forcing use of specific methods: The
178  // behavior (though not the speed) will be identical to invert(ifail).
179 
180  void invertHaywood4 (G4int & ifail);
181  void invertHaywood5 (G4int &ifail);
182  void invertHaywood6 (G4int &ifail);
183  void invertBunchKaufman (G4int &ifail);
184 
185  protected:
186 
187  inline G4int num_size() const;
188 
189  private:
190 
191  friend class G4ErrorSymMatrix_row;
193  friend class G4ErrorMatrix;
194 
195  friend void tridiagonal(G4ErrorSymMatrix *a, G4ErrorMatrix *hsm);
196  friend G4double condition(const G4ErrorSymMatrix &m);
197  friend void diag_step(G4ErrorSymMatrix *t, G4int begin, G4int end);
198  friend void diag_step(G4ErrorSymMatrix *t, G4ErrorMatrix *u,
199  G4int begin, G4int end);
202  G4int row, G4int col);
203 
204  friend G4ErrorSymMatrix operator+(const G4ErrorSymMatrix &m1,
205  const G4ErrorSymMatrix &m2);
206  friend G4ErrorSymMatrix operator-(const G4ErrorSymMatrix &m1,
207  const G4ErrorSymMatrix &m2);
208  friend G4ErrorMatrix operator*(const G4ErrorSymMatrix &m1,
209  const G4ErrorSymMatrix &m2);
210  friend G4ErrorMatrix operator*(const G4ErrorSymMatrix &m1,
211  const G4ErrorMatrix &m2);
212  friend G4ErrorMatrix operator*(const G4ErrorMatrix &m1,
213  const G4ErrorSymMatrix &m2);
214  // Multiply a Matrix by a Matrix or Vector.
215 
216  // Returns v * v.T();
217 
218  std::vector<G4double > m;
219 
221  G4int size; // total number of elements
222 
227 
232 
233  void invert4 (G4int & ifail);
234  void invert5 (G4int & ifail);
235  void invert6 (G4int & ifail);
236 };
237 
238 //
239 // Operations other than member functions for Matrix, G4ErrorSymMatrix,
240 // DiagMatrix and Vectors
241 //
242 
243 std::ostream& operator<<(std::ostream &s, const G4ErrorSymMatrix &q);
244 // Write out Matrix, G4ErrorSymMatrix, DiagMatrix and Vector into ostream.
245 
247  const G4ErrorSymMatrix &m2);
249  const G4ErrorMatrix &m2);
251  const G4ErrorSymMatrix &m2);
254 // Multiplication operators.
255 // Note that m *= m1 is always faster than m = m * m1
256 
258 // s = s1 / t. (s /= t is faster if you can use it.)
259 
261  const G4ErrorSymMatrix &s2);
263  const G4ErrorMatrix &m2);
265  const G4ErrorSymMatrix &s2);
266 // Addition operators
267 
269  const G4ErrorSymMatrix &s2);
271  const G4ErrorMatrix &m2);
273  const G4ErrorSymMatrix &s2);
274 // subtraction operators
275 
277  const G4ErrorSymMatrix &s2);
278 // Direct sum of two symmetric matrices;
279 
281 // Find the conditon number of a symmetric matrix.
282 
283 void diag_step(G4ErrorSymMatrix *t, G4int begin, G4int end);
284 void diag_step(G4ErrorSymMatrix *t, G4ErrorMatrix *u, G4int begin, G4int end);
285 // Implicit symmetric QR step with Wilkinson Shift
286 
288 // Diagonalize a symmetric matrix.
289 // It returns the matrix U so that s_old = U * s_diag * U.T()
290 
292  G4int row=1, G4int col=1);
293 // Finds and does Householder reflection on matrix.
294 
297 // Does a Householder tridiagonalization of a symmetric matrix.
298 
299 #include "G4ErrorSymMatrix.icc"
300 
301 #endif
void house_with_update2(G4ErrorSymMatrix *a, G4ErrorMatrix *v, G4int row=1, G4int col=1)
const G4double & fast(G4int row, G4int col) const
void invert4(G4int &ifail)
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
friend void tridiagonal(G4ErrorSymMatrix *a, G4ErrorMatrix *hsm)
static const G4double CHOLESKY_THRESHOLD_5x5
std::ostream & operator<<(std::ostream &, const BasicVector3D< float > &)
friend void diag_step(G4ErrorSymMatrix *t, G4int begin, G4int end)
void invertCholesky6(G4int &ifail)
void invertBunchKaufman(G4int &ifail)
G4ErrorSymMatrix_row(G4ErrorSymMatrix &, G4int)
static constexpr double m2
Definition: G4SIunits.hh:130
const G4double & operator()(G4int row, G4int col) const
friend G4ErrorMatrix operator*(const G4ErrorSymMatrix &m1, const G4ErrorSymMatrix &m2)
const char * p
Definition: xmltok.h:285
friend G4ErrorSymMatrix operator+(const G4ErrorSymMatrix &m1, const G4ErrorSymMatrix &m2)
G4ErrorSymMatrix apply(G4double(*f)(G4double, G4int, G4int)) const
virtual ~G4ErrorSymMatrix()
void invert(G4int &ifail)
friend G4double condition(const G4ErrorSymMatrix &m)
G4double condition(const G4ErrorSymMatrix &m)
void invert6(G4int &ifail)
#define G4ThreadLocal
Definition: tls.hh:69
static G4ThreadLocal G4double posDefFraction5x5
G4ErrorSymMatrix operator-() const
void invertCholesky5(G4int &ifail)
BasicVector3D< float > operator*(const BasicVector3D< float > &v, double a)
static G4ThreadLocal G4double posDefFraction6x6
static constexpr double m
Definition: G4SIunits.hh:129
const XML_Char * s
Definition: expat.h:262
double G4double
Definition: G4Types.hh:76
void tridiagonal(G4ErrorSymMatrix *a, G4ErrorMatrix *hsm)
friend G4ErrorMatrix diagonalize(G4ErrorSymMatrix *s)
G4ErrorMatrix dsum(const G4ErrorMatrix &, const G4ErrorMatrix &)
Int_t col[ntarg]
Definition: Style.C:29
BasicVector3D< float > operator+(const BasicVector3D< float > &v)
static G4ThreadLocal G4double adjustment6x6
G4double trace() const
G4ErrorSymMatrix similarity(const G4ErrorMatrix &m1) const
static G4ThreadLocal G4double adjustment5x5
G4ErrorSymMatrix & operator+=(const G4ErrorSymMatrix &m2)
const G4double & operator[](G4int) const
G4int num_row() const
std::vector< G4double > m
G4ErrorSymMatrix similarityT(const G4ErrorMatrix &m1) const
G4int num_col() const
friend void house_with_update2(G4ErrorSymMatrix *a, G4ErrorMatrix *v, G4int row, G4int col)
int G4int
Definition: G4Types.hh:78
void assign(const G4ErrorMatrix &m2)
static const G4double CHOLESKY_CREEP_6x6
G4ErrorSymMatrix & operator=(const G4ErrorSymMatrix &m2)
void invert5(G4int &ifail)
G4ErrorSymMatrix sub(G4int min_row, G4int max_row) const
G4ErrorSymMatrix_row_const(const G4ErrorSymMatrix &, G4int)
G4ErrorMatrix diagonalize(G4ErrorSymMatrix *s)
BasicVector3D< float > operator/(const BasicVector3D< float > &v, double a)
void invertHaywood4(G4int &ifail)
static const G4double CHOLESKY_THRESHOLD_6x6
G4double determinant() const
void diag_step(G4ErrorSymMatrix *t, G4int begin, G4int end)
static const G4double CHOLESKY_CREEP_5x5
G4ErrorSymMatrix T() const
void invertHaywood6(G4int &ifail)
void invertHaywood5(G4int &ifail)
G4int num_size() const
G4ErrorSymMatrix & operator-=(const G4ErrorSymMatrix &m2)
G4ErrorSymMatrix & operator*=(G4double t)
G4ErrorSymMatrix inverse(G4int &ifail) const
BasicVector3D< float > operator-(const BasicVector3D< float > &v)
G4ErrorSymMatrix_row operator[](G4int)
G4ErrorSymMatrix & operator/=(G4double t)