Jpp - the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JMatrix4S.cc
Go to the documentation of this file.
1 #include <string>
2 #include <iostream>
3 #include <iomanip>
4 
5 #include "TRandom3.h"
6 #include "TMatrixTSym.h"
7 
8 #include "JLang/JException.hh"
9 #include "JMath/JMatrix4S.hh"
10 #include "JMath/JMathTestkit.hh"
11 
12 #include "Jeep/JParser.hh"
13 #include "Jeep/JMessage.hh"
14 #include "Jeep/JTimer.hh"
15 
16 
17 /**
18  * \file
19  *
20  * Example program to test inversion of symmetric matrix.
21  * \author mdejong
22  */
23 int main(int argc, char**argv)
24 {
25  using namespace std;
26  using namespace JPP;
27 
28  double precision;
29  int numberOfEvents;
30  int debug;
31 
32  try {
33 
34  JParser<> zap("Example program to test inversion of symmetric matrix.");
35 
36  zap['e'] = make_field(precision) = 1.0e-6;
37  zap['n'] = make_field(numberOfEvents) = 1000;
38  zap['d'] = make_field(debug) = 3;
39 
40  zap(argc, argv);
41  }
42  catch(const exception &error) {
43  FATAL(error.what() << endl);
44  }
45 
46 
47  ASSERT(numberOfEvents > 0);
48 
49  gRandom->SetSeed(0);
50 
51  JTimer t1("Jpp");
52  JTimer t2("ROOT");
53 
54  for (int i = 0; i != numberOfEvents; ++i) {
55 
56  if (i%100 == 0) { STATUS(setw(8) << i << '\r'); DEBUG(endl); }
57 
58  JMatrix4S A = getRandom<JMatrix4S>();
59 
60  TMatrixTSym<double> R(4);
61 
62  R(0,0) = A.a00; R(0,1) = A.a01; R(0,2) = A.a02; R(0,3) = A.a03;
63  R(1,0) = A.a10; R(1,1) = A.a11; R(1,2) = A.a12; R(1,3) = A.a13;
64  R(2,0) = A.a20; R(2,1) = A.a21; R(2,2) = A.a22; R(2,3) = A.a23;
65  R(3,0) = A.a30; R(3,1) = A.a31; R(3,2) = A.a32; R(3,3) = A.a33;
66 
67 
68  DEBUG("Matrix A" << endl);
69  DEBUG(A << endl);
70 
71  DEBUG("Determinant A " << A.getDeterminant() << endl);
72 
73  try {
74 
75  JMatrix4S B(A);
76 
77  t1.start();
78 
79  for (int i = 11; i != 0; --i) {
80  B.invert();
81  }
82 
83  t1.stop();
84 
85  t2.start();
86 
87  for (int i = 11; i != 0; --i) {
88  R.Invert();
89  }
90 
91  t2.stop();
92 
93  DEBUG("Matrix A^-1" << endl);
94  DEBUG(B << endl);
95 
96  DEBUG("Determinant A^-1 = " << B.getDeterminant() << endl);
97 
98  JMatrix4D C(A * B);
99 
100  DEBUG("Matrix A x A^-1" << endl);
101  DEBUG(C << endl);
102 
103  DEBUG("Determinant (A x A^-1) = " << C.getDeterminant() << endl);
104  DEBUG("Determinant A x Determinant A^-1 = " << A.getDeterminant() * B.getDeterminant() << endl);
105  DEBUG("A x A^-1 = I ? " << C.isIdentity(precision) << endl);
106 
107  if (!C.isIdentity(precision)) {
108  ERROR("Matrix A x A^-1 /= I" << endl);
109  }
110 
111  JMatrix4D D = C - JMatrix4D::getIdentity();
112 
113  DEBUG("Matrix D = C - I" << endl);
114  DEBUG(D << endl);
115 
116  ASSERT(C.isIdentity(precision));
117  }
118  catch (const JException& error) {
119  FATAL(error << endl);
120  }
121  }
122  STATUS(endl);
123 
124  if (numberOfEvents > 0) {
125 
126  const double factor = 1.0; // / numberOfEvents;
127 
128  t1.print(cout, factor, micro_t);
129  t2.print(cout, factor, micro_t);
130  }
131 
132  return 0;
133 }
Utility class to parse command line options.
Definition: JParser.hh:1500
Exceptions.
do echo Generating $dir eval D
Definition: JDrawLED.sh:50
int main(int argc, char *argv[])
Definition: Main.cc:15
micro
Definition: JScale.hh:29
#define STATUS(A)
Definition: JMessage.hh:63
static const double C
Physics constants.
#define ASSERT(A,...)
Assert macro.
Definition: JMessage.hh:90
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1961
#define ERROR(A)
Definition: JMessage.hh:66
int debug
debug level
Definition: JSirene.cc:63
then usage $script[distance] fi case set_variable R
Definition: JDrawLED.sh:40
General purpose messaging.
#define FATAL(A)
Definition: JMessage.hh:67
Utility class to parse command line options.
source $JPP_DIR setenv csh $JPP_DIR eval JShellParser o a A
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62