CLHEP VERSION Reference Documentation
   
CLHEP Home Page     CLHEP Documentation     CLHEP Bug Reports

GenMatrix.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 // ---------------------------------------------------------------------------
3 //
4 // This file is a part of the CLHEP - a Class Library for High Energy Physics.
5 //
6 // This is the implementation of the HepGenMatrix class.
7 //
8 
9 #ifdef GNUPRAGMA
10 #pragma implementation
11 #endif
12 
13 #include <string.h>
14 #include <cmath>
15 #include <stdlib.h>
16 
17 #include "CLHEP/Matrix/GenMatrix.h"
18 #include "CLHEP/Matrix/SymMatrix.h"
19 #include "CLHEP/Matrix/Matrix.h"
20 
21 #ifdef HEP_DEBUG_INLINE
22 #include "CLHEP/Matrix/GenMatrix.icc"
23 #endif
24 
25 namespace CLHEP {
26 
27 #ifdef HEP_THIS_FUNCTION_IS_NOT_NEEDED
28 static void delete_array(double *hm)
29 {
30  delete [] hm;
31 }
32 #endif
33 
34 double norm_infinity(const HepGenMatrix &hm) {
35  double max=0,sum;
36  for(int r=1;r<=hm.num_row();r++) {
37  sum=0;
38  for(int c=1;c<=hm.num_col();c++) {
39  sum+=fabs(hm(r,c));
40  }
41  if(sum>max) max=sum;
42  }
43  return max;
44 }
45 
46 double norm1(const HepGenMatrix &hm) {
47  double max=0,sum;
48  for(int c=1;c<=hm.num_col();c++) {
49  sum=0;
50  for(int r=1;r<=hm.num_row();r++)
51  sum+=fabs(hm(r,c));
52  if(sum>max) max=sum;
53  }
54  return max;
55 }
56 
57 double norm(const HepGenMatrix &hm) {
58  HepSymMatrix A(hm.num_col(),0);
59 
60 // Calculate hm.T*hm
61  int r;
62  for(r=1;r<=A.num_row();r++)
63  for(int c=1;c<=r;c++)
64  for(int i=1;i<=hm.num_row();i++)
65  A.fast(r,c)=hm(i,r)*hm(i,c);
66  diagonalize(&A);
67  double max=fabs(A(1,1));
68  for(r=2;r<=A.num_row();r++)
69  if(max<fabs(A(r,r))) max=fabs(A(r,r));
70  return (sqrt(max));
71 }
72 
73 void HepGenMatrix::error(const char *es)
74 {
75  std::cerr << es << std::endl;
76  std::cerr << "---Exiting to System." << std::endl;
77  abort();
78 }
79 
80 bool HepGenMatrix::operator== ( const HepGenMatrix& o) const {
81  if(o.num_row()!=num_row() || o.num_col()!=num_col()) return false;
82  for (int k1=1; k1<=num_row(); k1++)
83  for (int k2=1; k2<=num_col(); k2++)
84  if(o(k1,k2) != (*this)(k1,k2)) return false;
85  return true;
86 }
87 
88 // implementation using pre-allocated data array
89 // -----------------------------------------------------------------
90 
91 void HepGenMatrix::delete_m(int size, double* hm)
92 {
93  if (hm)
94  {
95  if(size > size_max)
96  delete [] hm;
97  }
98 }
99 
100 double* HepGenMatrix::new_m(int )
101 {
102  /*-ap: data_array is replaced by the std::vector<double>,
103  * so we simply return 0 here
104  *
105  * if (size == 0) return 0;
106  * else {
107  * if ( size <= size_max ) {
108  * memset(data_array, 0, size * sizeof(double));
109  * return data_array;
110  * } else {
111  * double * nnn = new double[size];
112  * memset(nnn, 0, size * sizeof(double));
113  * return nnn;
114  * }
115  * }
116  *-ap end
117  */
118  return 0;
119 }
120 
121 } // namespace CLHEP
122 
double norm(const HepGenMatrix &m)
Definition: GenMatrix.cc:57
virtual int num_row() const =0
#define abort()
double norm1(const HepGenMatrix &m)
Definition: GenMatrix.cc:46
virtual int num_col() const =0
virtual bool operator==(const HepGenMatrix &) const
Definition: GenMatrix.cc:80
const double & fast(int row, int col) const
void delete_m(int size, double *)
Definition: GenMatrix.cc:91
double * new_m(int size)
Definition: GenMatrix.cc:100
static void error(const char *s)
Definition: GenMatrix.cc:73
double norm_infinity(const HepGenMatrix &m)
Definition: GenMatrix.cc:34
HepMatrix diagonalize(HepSymMatrix *s)