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