ATLAS Offline Software
Public Member Functions | Protected Member Functions | Static Protected Member Functions | Protected Attributes | List of all members
Trk::AlSpaMat Class Reference

#include <AlSpaMat.h>

Inheritance diagram for Trk::AlSpaMat:
Collaboration diagram for Trk::AlSpaMat:

Public Member Functions

 AlSpaMat ()
 
 AlSpaMat (long int N)
 
 AlSpaMat (const AlSpaMat &m)
 
 AlSpaMat (const AlSymMat &m)
 
 ~AlSpaMat ()
 
AlSpaMatoperator= (const double &)
 
AlSpaMatoperator= (const AlSpaMat &)
 
AlSpaMatoperator= (const AlSymMat &)
 
AlSpaMatoperator= (const AlMat &)
 
AlSpaMat operator+ (const AlSpaMat &) const
 
AlSpaMatoperator+= (const AlSpaMat &)
 
AlSpaMat operator- (const AlSpaMat &) const
 
AlSpaMatoperator-= (const AlSpaMat &)
 
AlMat operator* (const AlSymMatBase &) const
 
AlVec operator* (const AlVec &) const
 
AlMat operator* (const AlMat &) const
 
AlSpaMat operator* (const double &) const
 
AlSpaMatoperator*= (const double &)
 
int SolveWithEigen (AlVec &RHS)
 
virtual void RemoveModule (int) override final
 
virtual void RemoveAlignPar (int, int) override final
 
void RemoveDoF (int, int nelem=1)
 
virtual int RemoveCollsRows (std::vector< int >) override final
 
virtual void reSize (long int) override final
 
virtual void SetPathBin (const std::string &) override final
 
virtual void SetPathTxt (const std::string &) override final
 
virtual StatusCode Write (const std::string &, bool, bool, double, float) override final
 
virtual StatusCode CheckMatVersion (const std::string &, bool &) override final
 
virtual StatusCode Read (const std::string &, int &, bool &, float &) override final
 
virtual StatusCode ReadProjected (const std::string &, int &, bool &, float &) override final
 
virtual int invert () override final
 
virtual int diagonalize (char jobz, AlVec &w, AlMat &z) override final
 
virtual double & elemr (long int, long int) override final
 
virtual double elemc (long int, long int) const override final
 
long int nele ()
 
virtual double determinant () override final
 
void elem (const indices &, long int &, long int &) const
 
const std::string & pathBin () const
 
const std::string & pathTxt () const
 
virtual TMatrixDSparsemakeTMatrix () override final
 
AlSymMatBase_row operator[] (long int)
 
AlSymMatBase_row_const operator[] (long int) const
 
long int nrow () const
 
long int ncol () const
 
long int size () const
 
int matrix_type () const
 
const datamapptrMap () const
 

Protected Member Functions

void copy (const AlSpaMat &m)
 
void copy (const AlSymMat &m)
 
void copy (const AlMat &m)
 

Static Protected Member Functions

static indices elem (long int, long int)
 

Protected Attributes

std::string m_pathbin
 
std::string m_pathtxt
 
int * m_ptr_row = nullptr
 
int * m_ptr_col = nullptr
 
int m_matrix_type = 0
 
datamap m_ptr_map
 
long int m_size = 0
 
long int m_nele = 0
 

Detailed Description

contains the implementation for handling sparse matrices

Definition at line 27 of file AlSpaMat.h.

Constructor & Destructor Documentation

◆ AlSpaMat() [1/4]

Trk::AlSpaMat::AlSpaMat ( )

Definition at line 43 of file AlSpaMat.cxx.

44  : m_ptr_row(nullptr)
45  , m_ptr_col(nullptr)
46 {
47  m_matrix_type = 2;
48  m_size = 0;
49  m_nele = 0;
50 }

◆ AlSpaMat() [2/4]

Trk::AlSpaMat::AlSpaMat ( long int  N)

Definition at line 53 of file AlSpaMat.cxx.

54  : m_ptr_row(nullptr)
55  , m_ptr_col(nullptr)
56 {
57  m_matrix_type = 2;
58  m_size = N;
59  m_nele = 0;
60 }

◆ AlSpaMat() [3/4]

Trk::AlSpaMat::AlSpaMat ( const AlSpaMat m)

Definition at line 63 of file AlSpaMat.cxx.

64  : AlSymMatBase(m)
65  , m_ptr_row(nullptr)
66  , m_ptr_col(nullptr)
67 {
68  m_matrix_type = 2;
69  m_size = m.size();
70  copy(m);
71 }

◆ AlSpaMat() [4/4]

Trk::AlSpaMat::AlSpaMat ( const AlSymMat m)

Definition at line 74 of file AlSpaMat.cxx.

75  : m_ptr_row(nullptr)
76  , m_ptr_col(nullptr)
77 {
78  m_matrix_type = 2;
79  m_size = m.size();
80  copy(m);
81 }

◆ ~AlSpaMat()

Trk::AlSpaMat::~AlSpaMat ( )

Definition at line 84 of file AlSpaMat.cxx.

85 {
86  m_ptr_map.clear();
87  if( m_ptr_row != nullptr ) delete [] m_ptr_row;
88  if( m_ptr_col != nullptr ) delete [] m_ptr_col;
89 }

Member Function Documentation

◆ CheckMatVersion()

StatusCode Trk::AlSpaMat::CheckMatVersion ( const std::string &  filename,
bool &  StdUnits 
)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 712 of file AlSpaMat.cxx.

713 {
714  std::ifstream inmat((filename).c_str(), std::ios::binary);
715  if(inmat.fail())
716  return StatusCode::FAILURE;
717 
718  m_ptr_map.clear();
719 
720  int32_t msiz=0;
721  inmat.read(charAddress(msiz), sizeof (msiz));
722 
723  float version=0.0;
724  inmat.read(charAddress(version), sizeof (version));
725 
726  StdUnits = version>=2.0;
727 
728  inmat.close();
729 
730  return StatusCode::SUCCESS;
731 }

◆ copy() [1/3]

void Trk::AlSpaMat::copy ( const AlMat m)
protected

Definition at line 122 of file AlSpaMat.cxx.

123 {
124  if( size() != m.nrow() || size() != m.ncol() ) {
125  throw std::range_error( "AlSpaMat::copy: size do not match!" );
126  }
127 
128  // copy just the lower triangle:
129  m_ptr_map.clear();
130  m_nele=0;
131 
132  //Convert Storage System
133  for (int i = 0; i<m_size; i++)
134  for (int j = 0; j<=i; j++)
135  if (m.elemc(i,j) != 0.) {
136  m_ptr_map.insert(std::make_pair(std::make_pair(i,j), m.elemc(i,j)));
137  m_nele++;
138  }
139 
140  }

◆ copy() [2/3]

void Trk::AlSpaMat::copy ( const AlSpaMat m)
protected

Definition at line 92 of file AlSpaMat.cxx.

93 {
94  if( size() != m.size()) {
95  throw std::range_error( "AlSpaMat::copy: size do not match!" );
96  }
97  m_ptr_map.clear();
98  m_nele=m.m_nele;
99  m_ptr_map = m.m_ptr_map;
100 }

◆ copy() [3/3]

void Trk::AlSpaMat::copy ( const AlSymMat m)
protected

Definition at line 103 of file AlSpaMat.cxx.

104 {
105  if( size() != m.size()) {
106  throw std::range_error( "AlSpaMat::copy: size do not match!" );
107  }
108  m_ptr_map.clear();
109  m_nele=0;
110 
111  // Convert Storage System
112  for (int i = 0; i<m_size; i++)
113  for (int j = 0; j<=i; j++)
114  if (m.elemc(i,j) != 0.) {
115  m_ptr_map.insert(std::make_pair(std::make_pair(i,j), m.elemc(i,j)));
116  m_nele++;
117  }
118 
119  }

◆ determinant()

double Trk::AlSpaMat::determinant ( )
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 403 of file AlSpaMat.cxx.

404 {
405  throw std::invalid_argument( "AlSpaMat::determinant: not implemented!" );
406 }

◆ diagonalize()

int Trk::AlSpaMat::diagonalize ( char  jobz,
AlVec w,
AlMat z 
)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 481 of file AlSpaMat.cxx.

482 {
483  throw std::invalid_argument( "AlSpaMat::diagonalize: not implemented!" );
484 }

◆ elem() [1/2]

void Trk::AlSpaMat::elem ( const indices key,
long int &  i,
long int &  j 
) const
inline

Definition at line 94 of file AlSpaMat.h.

94  {
95  i = key.first;
96  j = key.second;
97  if(j>i) throw std::out_of_range( "AlSpaMat::elem: i<j ! " );
98  return;
99 }

◆ elem() [2/2]

indices Trk::AlSpaMat::elem ( long int  i,
long int  j 
)
staticprotected

Definition at line 209 of file AlSpaMat.cxx.

210 {
211  // ATTENTION! the key value is returned:
212 #ifdef _DEBUG
213  if( i<0 ) {
214  throw std::out_of_range( "AlSpaMat::elem: Index 1 < zero! " );
215  }
216  if( i>=size() ) {
217  throw std::out_of_range( "AlSpaMat::elem: Index 1 too large! " );
218  }
219  if( j<0 ) {
220  throw std::out_of_range( "AlSpaMat::elem: Index 2 < zero! " );
221  }
222  if( j>=size() ) {
223  throw std::out_of_range( "AlSpaMat::elem: Index 2 too large! " );
224  }
225 #endif
226 
227  indices key = j < i ? std::make_pair(i,j) : std::make_pair(j,i);
228  return key;
229 }

◆ elemc()

double Trk::AlSpaMat::elemc ( long int  i,
long int  j 
) const
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 179 of file AlSpaMat.cxx.

180 {
181 #ifdef _DEBUG
182  if( i<0 ) {
183  throw std::out_of_range( "AlSpaMat::elemc: Index 1 < zero! " );
184  }
185  if( i>=size() ) {
186  throw std::out_of_range( "AlSpaMat::elemc: Index 1 too large! " );
187  }
188  if( j<0 ) {
189  throw std::out_of_range( "AlSpaMat::elemc: Index 2 < zero! " );
190  }
191  if( j>=size() ) {
192  throw std::out_of_range( "AlSpaMat::elemc: Index 2 too large! " );
193  }
194 #endif
195  // try fast referencing:
197  indices key = j < i ? std::make_pair(i,j) : std::make_pair(j,i);
198  pos=m_ptr_map.find(key);
199 
200  if( pos!=m_ptr_map.end() ) // already exists
201  return pos->second;
202 
203  // does not yet exist
204  // return zero:
205  return 0.;
206 }

◆ elemr()

double & Trk::AlSpaMat::elemr ( long int  i,
long int  j 
)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 143 of file AlSpaMat.cxx.

144 {
145 #ifdef _DEBUG
146  if( i<0 ) {
147  throw std::out_of_range( "AlSpaMat::elemr: Index 1 < zero! " );
148  }
149  if( i>=size() ) {
150  throw std::out_of_range( "AlSpaMat::elemr: Index 1 too large! " );
151  }
152  if( j<0 ) {
153  throw std::out_of_range( "AlSpaMat::elemr: Index 2 < zero! " );
154  }
155  if( j>=size() ) {
156  throw std::out_of_range("AlSpaMat::elemr: Index 2 too large! " );;
157  }
158 #endif
159  // try fast referencing:
161  indices key = j < i ? std::make_pair(i,j) : std::make_pair(j,i);
162  pos=m_ptr_map.find(key);
163  if( pos!=m_ptr_map.end() ) // already exists
164  return pos->second;
165  else { // does not yet exist
166  // create it with the null value:
167  m_nele++;
168  return (m_ptr_map.insert(std::make_pair(key, 0.))).first->second;
169  }
170  // simpler implementation of the above:
171  /*
172  m_ptr_map[key]=0;
173  m_nele = m_ptr_map.size();
174  return m_ptr_map[key];
175  */
176 }

◆ invert()

int Trk::AlSpaMat::invert ( )
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 487 of file AlSpaMat.cxx.

488 {
489  throw std::invalid_argument( "AlSpaMat::invert: not implemented!" );
490 }

◆ makeTMatrix()

TMatrixDSparse * Trk::AlSpaMat::makeTMatrix ( )
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 861 of file AlSpaMat.cxx.

862 {
863 
864  long int nonZeroElements = m_ptr_map.size();
865  int *irow = new int[nonZeroElements];
866  int *icol = new int[nonZeroElements];
867  double *val = new double[nonZeroElements];
868 
869  long int i, j;
870  long int counter(0);
871  for (const datamap::value_type& p : m_ptr_map) {
872  i = p.first.first;
873  j = p.first.second;
874  *(val+counter)=p.second;
875  *(irow+counter)= i;
876  *(icol+counter)= j;
877  counter++;
878  }
879 
880  TMatrixDSparse* myTMatrix = new TMatrixDSparse(0,m_size-1,0,m_size-1);
881  myTMatrix->SetMatrixArray(nonZeroElements,irow,icol,val);
882 
883  delete [] irow;
884  delete [] icol;
885  delete [] val;
886  return myTMatrix;
887 }

◆ matrix_type()

int Trk::AlSymMatBase::matrix_type ( ) const
inlineinherited

Definition at line 148 of file AlSymMatBase.h.

148 { return m_matrix_type; }

◆ ncol()

long int Trk::AlSymMatBase::ncol ( ) const
inlineinherited

Definition at line 147 of file AlSymMatBase.h.

147 { return m_size; }

◆ nele()

long int Trk::AlSpaMat::nele ( )
inline

Definition at line 101 of file AlSpaMat.h.

101  {
102  // first check for intrinsic consistency:
103  if( int(m_ptr_map.size()) != m_nele ) throw std::range_error( "AlSpaMat::m_nele has been corrupted!" );
104  m_nele = m_ptr_map.size();
105  return m_nele;
106 }

◆ nrow()

long int Trk::AlSymMatBase::nrow ( ) const
inlineinherited

Definition at line 146 of file AlSymMatBase.h.

146 { return m_size; }

◆ operator*() [1/4]

AlMat Trk::AlSpaMat::operator* ( const AlMat m) const

Definition at line 349 of file AlSpaMat.cxx.

349  {
350  if( size() != m.nrow() ) {
351  throw std::range_error( "AlSpaMat::operator*: size do not match!" );
352  }
353 
354  long int isiz(size());
355  long int icol(m.ncol());
356  AlMat b( isiz, icol);
357 
358  for( long int i=0; i<isiz; i++ )
359  for( long int j=0; j<icol; j++ )
360  for( long int k=0; k<isiz; k++)
361  b.elemr(i,j) += elemc(i,k)*m.elemc(k,j);
362 
363  return b;
364 }

◆ operator*() [2/4]

AlMat Trk::AlSpaMat::operator* ( const AlSymMatBase m) const

Definition at line 331 of file AlSpaMat.cxx.

332 {
333  if( size() != m.size() ) {
334  throw std::range_error( "AlSpaMat::operator*: size do not match!" );
335  }
336 
337  long int isiz(size());
338  AlMat b( isiz, isiz);
339 
340  for( long int i=0; i<isiz; i++ )
341  for( long int j=0; j<isiz; j++ )
342  for( int k=0; k<isiz; k++)
343  b.elemr(i,j) += elemc(i,k)*m.elemc(k,j);
344 
345  return b;
346 }

◆ operator*() [3/4]

AlVec Trk::AlSpaMat::operator* ( const AlVec v) const

Definition at line 367 of file AlSpaMat.cxx.

368 {
369  if( size() != v.size() ) {
370  throw std::range_error( "AlSpaMat::operator*: size do not match! " );
371  }
372 
373  long int isiz(size());
374  AlVec b(isiz);
375 
376  for( long int i=0; i<isiz; i++ )
377  for( long int j=0; j<isiz; j++ )
378  b[i] += elemc(i,j)*v[j];
379 
380  return b;
381 }

◆ operator*() [4/4]

AlSpaMat Trk::AlSpaMat::operator* ( const double &  d) const

Definition at line 393 of file AlSpaMat.cxx.

394 {
395  AlSpaMat a(size());
396  for (const datamap::value_type& p : m_ptr_map)
397  a.m_ptr_map.emplace(p.first, p.second*d);
398 
399  return a;
400 }

◆ operator*=()

AlSpaMat & Trk::AlSpaMat::operator*= ( const double &  d)

Definition at line 384 of file AlSpaMat.cxx.

385 {
386  for (datamap::value_type& p : m_ptr_map)
387  p.second *= d;
388 
389  return *this;
390 }

◆ operator+()

AlSpaMat Trk::AlSpaMat::operator+ ( const AlSpaMat m) const

Definition at line 273 of file AlSpaMat.cxx.

274 {
275  if( size() != m.size()) {
276  throw std::range_error( "AlSpaMat::operator+: size do not match!" );
277  }
278 
279  AlSpaMat b(m);
280  for (const datamap::value_type& p : m_ptr_map)
281  b.m_ptr_map[p.first] += p.second;
282  b.m_nele = b.m_ptr_map.size();
283 
284  return b;
285 }

◆ operator+=()

AlSpaMat & Trk::AlSpaMat::operator+= ( const AlSpaMat m)

Definition at line 288 of file AlSpaMat.cxx.

289 {
290  if( size() != m.size()) {
291  throw std::range_error( "AlSpaMat::operator+=: size do not match!" );
292  }
293 
294  for (const datamap::value_type& p : m.m_ptr_map)
295  (*this).m_ptr_map[p.first] += p.second;
296  m_nele = m_ptr_map.size();
297 
298  return *this;
299 }

◆ operator-()

AlSpaMat Trk::AlSpaMat::operator- ( const AlSpaMat m) const

Definition at line 302 of file AlSpaMat.cxx.

303 {
304  if( size() != m.size()) {
305  throw std::range_error( "AlSpaMat::operator-: size do not match!" );
306  }
307 
308  AlSpaMat b(m);
309  for (const datamap::value_type& p : m_ptr_map)
310  b.m_ptr_map[p.first] -= p.second;
311  b.m_nele = b.m_ptr_map.size();
312 
313  return b;
314 }

◆ operator-=()

AlSpaMat & Trk::AlSpaMat::operator-= ( const AlSpaMat m)

Definition at line 317 of file AlSpaMat.cxx.

318 {
319  if( size() != m.size()) {
320  throw std::range_error( "AlSpaMat::operator-=: size do not match!" );
321  }
322 
323  for (const datamap::value_type& p : m.m_ptr_map)
324  (*this).m_ptr_map[p.first] -= p.second;
325  m_nele = m_ptr_map.size();
326 
327  return *this;
328 }

◆ operator=() [1/4]

AlSpaMat & Trk::AlSpaMat::operator= ( const AlMat m)

Definition at line 252 of file AlSpaMat.cxx.

253 {
254  if( m.nrow() != m.ncol() ) {
255  throw std::range_error( "AlSpaMat::=operator: allowed for square matrices only!" );
256  }
257 
258  m_size=m.nrow();
259  copy(m);
260  return *this;
261 }

◆ operator=() [2/4]

AlSpaMat & Trk::AlSpaMat::operator= ( const AlSpaMat m)

Definition at line 233 of file AlSpaMat.cxx.

234 {
235 
236  if ( this != &m ) {
237  m_size=m.size();
238  copy(m);
239  }
240  return *this;
241 }

◆ operator=() [3/4]

AlSpaMat & Trk::AlSpaMat::operator= ( const AlSymMat m)

Definition at line 244 of file AlSpaMat.cxx.

245 {
246  m_size=m.size();
247  copy(m);
248  return *this;
249 }

◆ operator=() [4/4]

AlSpaMat & Trk::AlSpaMat::operator= ( const double &  d)

Definition at line 264 of file AlSpaMat.cxx.

265 {
266  for (datamap::value_type& p : m_ptr_map)
267  p.second = d;
268 
269  return *this;
270 }

◆ operator[]() [1/2]

AlSymMatBase::AlSymMatBase_row Trk::AlSymMatBase::operator[] ( long int  r)
inlineinherited

Definition at line 109 of file AlSymMatBase.h.

109  {
110  AlSymMatBase_row b(*this,r);
111  return b;
112 }

◆ operator[]() [2/2]

AlSymMatBase::AlSymMatBase_row_const Trk::AlSymMatBase::operator[] ( long int  r) const
inlineinherited

Definition at line 114 of file AlSymMatBase.h.

114  {
115  const AlSymMatBase_row_const b(*this,r);
116  return b;
117 }

◆ pathBin()

const std::string & Trk::AlSpaMat::pathBin ( ) const
inline

Definition at line 108 of file AlSpaMat.h.

108  {
109  return m_pathbin;
110 }

◆ pathTxt()

const std::string & Trk::AlSpaMat::pathTxt ( ) const
inline

Definition at line 112 of file AlSpaMat.h.

112  {
113  return m_pathtxt;
114 }

◆ ptrMap()

const datamap * Trk::AlSymMatBase::ptrMap ( ) const
inlineinherited

Definition at line 149 of file AlSymMatBase.h.

149 { return &m_ptr_map; }

◆ Read()

StatusCode Trk::AlSpaMat::Read ( const std::string &  filename,
int &  dofs,
bool &  triang,
float &  version 
)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 734 of file AlSpaMat.cxx.

735 {
736  bool stdUnits = true;
737  if (StatusCode::SUCCESS != CheckMatVersion(m_pathbin+filename, stdUnits))
738  return StatusCode::FAILURE;
739 
740  std::ifstream inmat((m_pathbin+filename).c_str(), std::ios::binary);
741  if(inmat.fail())
742  return StatusCode::FAILURE;
743 
744  m_ptr_map.clear();
745 
746  int32_t msiz=0;
747  int32_t nelem;
748  inmat.read(charAddress(msiz), sizeof (msiz));
749  if( msiz>999999 )
750  dofs = msiz-1000000;
751  else
752  dofs = abs(msiz);
753  m_size = dofs;
754 
755  if (stdUnits)
756  inmat.read(charAddress(version), sizeof (version));
757 
758  double melem=0;
759  int32_t i, j;
760 
761  if(msiz>999999) { // sparse format
762  triang=false;
763  inmat.read(charAddress(nelem), sizeof (nelem));
764  m_nele=nelem;
765  for(int k=0; k<nelem; k++) {
766  inmat.read(charAddress(i), sizeof (i));
767  inmat.read(charAddress(j), sizeof (j));
768  inmat.read(charAddress(melem), sizeof (melem));
769  m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
770  }
771  }
772  else if(msiz>0) { // square format
773  triang=false;
774  for(int32_t i=0; i<msiz; i++) {
775  for(int32_t j=0; j<msiz; j++) {
776  inmat.read(charAddress(melem), sizeof (melem));
777  if( i>=j && melem!=0. )
778  m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
779  }
780  }
781  }
782  else { // triangular format
783  triang=true;
784  msiz = (-1)*msiz;
785  for( int32_t i=0; i<msiz; i++) {
786  for( int32_t j=0; j<=i; j++) {
787  inmat.read(charAddress(melem), sizeof (melem));
788  if( melem!=0. )
789  m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
790  }
791  }
792  }
793 
794  m_nele = m_ptr_map.size();
795  inmat.close();
796  return StatusCode::SUCCESS;
797 }

◆ ReadProjected()

StatusCode Trk::AlSpaMat::ReadProjected ( const std::string &  filename,
int &  dofs,
bool &  triang,
float &  version 
)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 800 of file AlSpaMat.cxx.

802 {
803  std::ifstream inmat((m_pathbin+filename).c_str(), std::ios::binary);
804  if(inmat.fail())
805  return StatusCode::FAILURE;
806 
807  m_ptr_map.clear();
808 
809  int32_t msiz=0;
810  int32_t nelem;
811  inmat.read(charAddress(msiz), sizeof (msiz));
812  if( msiz>999999 )
813  dofs = msiz-1000000;
814  else
815  dofs = std::abs(msiz);
816  m_size = dofs;
817 
818  inmat.read(charAddress(version), sizeof (version));
819 
820  double melem=0;
821  int32_t i, j;
822 
823  if(msiz>999999) { // sparse format
824  triang=false;
825  inmat.read(charAddress(nelem), sizeof (nelem));
826  m_nele=nelem;
827  for(int k=0; k<nelem; k++) {
828  inmat.read(charAddress(i), sizeof (i));
829  inmat.read(charAddress(j), sizeof (j));
830  inmat.read(charAddress(melem), sizeof (melem));
831  m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
832  }
833  }
834  else if(msiz>0) { // square format
835  triang=false;
836  for(int32_t i=0; i<msiz; i++) {
837  for(int32_t j=0; j<msiz; j++) {
838  inmat.read(charAddress(melem), sizeof (melem));
839  if( i>=j && melem!=0. )
840  m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
841  }
842  }
843  }
844  else { // triangular format
845  triang=true;
846  msiz = (-1)*msiz;
847  for( int32_t i=0; i<msiz; i++) {
848  for( int32_t j=0; j<=i; j++) {
849  inmat.read(charAddress(melem), sizeof (melem));
850  if( melem!=0. )
851  m_ptr_map.insert(std::make_pair(std::make_pair(i,j), melem));
852  }
853  }
854  }
855 
856  m_nele = m_ptr_map.size();
857  inmat.close();
858  return StatusCode::SUCCESS;
859 }

◆ RemoveAlignPar()

void Trk::AlSpaMat::RemoveAlignPar ( int  index,
int  control 
)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 564 of file AlSpaMat.cxx.

565 {
566  // Dynamic shiftRow
567  int shiftRow = 1;
568  int counterRow = 0;
569 
570  // Dynamic shiftCol
571  int shiftCol;
572  int counterCol = 0;
573 
574  index = index-control;
575 
577 
578  for(int row=index; row<(size()-shiftRow); row++) {
579  shiftCol = 0;
580  counterCol = 0;
581 
582  for(int col=0; col<row+1; col++) {
583 
584  if (col==index) {
585  shiftCol = 1;
586  counterCol=0;
587  }
588 
589  pos = m_ptr_map.find(std::make_pair(row+shiftRow,col+shiftCol));
590  if( pos!=m_ptr_map.end() ) { // exists
591  elemr(row,col) = pos->second;
592  m_ptr_map.erase(pos);
593  }
594  else {
595  pos = m_ptr_map.find(std::make_pair(row,col));
596  if( pos!=m_ptr_map.end() )
597  m_ptr_map.erase(pos);
598  }
599  counterCol++;
600  if (counterCol==5-control) {
601  counterCol=0;
602  shiftCol++;
603  }
604  }
605 
606  counterRow++;
607  if (counterRow==5-control) {
608  counterRow=0;
609  shiftRow++;
610  }
611  }
612 
613  m_nele = m_ptr_map.size();
614 }

◆ RemoveCollsRows()

int Trk::AlSpaMat::RemoveCollsRows ( std::vector< int >  indices)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 531 of file AlSpaMat.cxx.

532 {
533  int n = indices.size();
534  if (n==0) {
535  return m_size;
536  }
537  if (n>m_size) {
538  throw std::invalid_argument( "AlSpaMat::RemoveCollsRows: Vector of indices larger than matrix size." );
539  return m_size;
540  }
541 
542  // first sort the list of indices descending
543  // maybe not the fastest way to do that but it works
544  for (int i=0;i<n-1;i++)
545  for (int j=i+1; j<n; j++)
546  if (indices[j] > indices[i]) {
547  int itmp = indices[i];
548  indices[i] = indices[j];
549  indices[j] = itmp;
550  }
551 
552  // remove rows and columns starting from largest indices
553  for (int i=0;i<n;i++) {
554  if (indices[i] > m_size-1) {
555  throw std::out_of_range( "AlSpaMat::RemoveCollsRows: Index goes beyond matrix." );
556  }
557  RemoveDoF(indices[i]);
558  }
559 
560  return m_size;
561 }

◆ RemoveDoF()

void Trk::AlSpaMat::RemoveDoF ( int  index,
int  nelem = 1 
)

Definition at line 493 of file AlSpaMat.cxx.

494 {
495  // index here is the first alignment parameter to remove, nelem number of consecutive ones
496 
497  m_size-=nelem;
498  long int n=index+nelem-1;
500  mapiterator pos_obs=m_ptr_map.end();
501  for (pos = m_ptr_map.begin(); pos!=m_ptr_map.end(); ++pos) {
502 
503  if( pos_obs!=m_ptr_map.end() )
504  m_ptr_map.erase(pos_obs);
505 
506  pos_obs = m_ptr_map.end();
507 
508  long int k, l;
509  elem(pos->first, k, l);
510 
511  if( (k>=index && k<=n) || (l>=index && l<=n) )
512  pos_obs = pos;
513  else {
514  if( k > n )
515  k -= nelem;
516  if( l > n )
517  l -= nelem;
518  if( k >= index || l >= index ) {
519  elemr(k, l) = pos->second;
520  pos_obs = pos;
521  }
522  }
523  }
524 
525  if( pos_obs!=m_ptr_map.end() )
526  m_ptr_map.erase(pos_obs);
527  m_nele = m_ptr_map.size();
528 }

◆ RemoveModule()

void Trk::AlSpaMat::RemoveModule ( int  index)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 617 of file AlSpaMat.cxx.

618 {
619  // index here is module number within matrix (size=6*max_index)
620  const int shift=6;
621  const int ind=shift*index;
622 
623  RemoveDoF(ind, shift);
624 }

◆ reSize()

void Trk::AlSpaMat::reSize ( long int  n)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 627 of file AlSpaMat.cxx.

628 {
629  // balanced binary tree!
630  if( n!=m_size ) {
631  // make a temporary copy:
632  AlSpaMat m(*this);
633  m_ptr_map.clear();
634  m_size = n;
635  long int i, j;
636  for (const datamap::value_type& p : m.m_ptr_map) {
637  m.elem(p.first, i, j);
638  if( i<n && j<n )
639  m_ptr_map.insert(p);
640  }
641  }
642 
643  m_nele = m_ptr_map.size();
644 }

◆ SetPathBin()

void Trk::AlSpaMat::SetPathBin ( const std::string &  path)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 647 of file AlSpaMat.cxx.

648 {
649  m_pathbin.assign(path);
650 }

◆ SetPathTxt()

void Trk::AlSpaMat::SetPathTxt ( const std::string &  path)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 653 of file AlSpaMat.cxx.

654 {
655  m_pathtxt.assign(path);
656 }

◆ size()

long int Trk::AlSymMatBase::size ( ) const
inlineinherited

Definition at line 145 of file AlSymMatBase.h.

145 { return m_size; }

◆ SolveWithEigen()

int Trk::AlSpaMat::SolveWithEigen ( AlVec RHS)

Definition at line 408 of file AlSpaMat.cxx.

408  {
409 
410  if(RHS.size() != size() ){
411  throw std::range_error( "AlSpaMat::SolveWithEigen vector size is incorrect" );
412  }
413 
414  Eigen::VectorXd eigenBigVector( RHS.size() );
415  for(int i(0); i< RHS.size(); ++i ){
416  eigenBigVector[i] = RHS[i];
417  }
418 
419  Eigen::SparseMatrix<double> eigenBigMatrix( m_size, m_size );
420 
421  using Triplet = Eigen::Triplet<double>;
422  std::vector<Triplet> tripletList;
423  tripletList.reserve(m_nele);
424  long int i, j;
425  for (const datamap::value_type& p : m_ptr_map) {
426  elem(p.first, i, j);
427  tripletList.emplace_back(i,j,p.second);
428  if(i!=j) tripletList.emplace_back(j,i,p.second);
429  }
430  eigenBigMatrix.setFromTriplets(tripletList.begin(), tripletList.end());
431 
432  // Eigen::CholmodSupernodalLLT is much much quicker (x50) so it would be great to move to that in the future
433  // requires an external package SuiteSparse
434  // BiCGSTAB was the fastest iterative solver in Eigen from a quick test
435  // SimplicialLDLT was the fastest direct solver in Eigen (x2 slower than BiCGSTAB )
436 
437  // Eigen::BiCGSTAB<Eigen::SparseMatrix<double> > solver;
438  Eigen::SimplicialLDLT<Eigen::SparseMatrix<double> > solver;
439 
440  solver.compute(eigenBigMatrix);
441  if(solver.info()!=Eigen::Success) {
442  // decomposition failed
443  throw std::domain_error("AlSpaMat::SolveWithEigen: failed to compute: -- is the input matrix singular?" );
444  return 1;
445  }
446 
447  Eigen::VectorXd x = solver.solve( eigenBigVector );
448  if(solver.info()!=Eigen::Success) {
449  // solving failed
450  throw std::domain_error("AlSpaMat::SolveWithEigen: Failed to solve: -- is your matrix singular? ");
451  return 2;
452  }
453 
454  //Copy results into vector
455  for(int i(0); i< RHS.size(); ++i ){
456  RHS[i] = x[i];
457  }
458 
459  //Check accuracy
460  Eigen::VectorXd residual = eigenBigMatrix * x - eigenBigVector;
461  double sumresidual = 0;
462  for( int i=0; i<residual.size(); ++i){
463  sumresidual += fabs(residual[i]);
464  }
465  sumresidual /= (double) residual.size();
466 
467  if( sumresidual > 1e-3 ){
468  throw std::overflow_error( "AlSpaMat::SolveWithEigen: your solution is no good! ");
469  return 3;
470  }
471 
472  return 0;
473 
474 
475 }

◆ Write()

StatusCode Trk::AlSpaMat::Write ( const std::string &  filename,
bool  binary,
bool  ,
double  ,
float  version 
)
finaloverridevirtual

Implements Trk::AlSymMatBase.

Definition at line 661 of file AlSpaMat.cxx.

663 {
664  std::ofstream outmat;
665 
666  int32_t msizz = 1000000+size();
667  int32_t nelem = m_ptr_map.size();
668 
669  if(binary) {
670  outmat.open((m_pathbin+filename).c_str(), std::ios::binary);
671  if(outmat.fail())
672  return StatusCode::FAILURE;
673  // add 10^6 to the size to distinguish
674  // from dense format
675  outmat.write(charAddress(msizz), sizeof (msizz));
676  outmat.write(charAddress(version), sizeof (version));
677  outmat.write(charAddress(nelem), sizeof (nelem));
678  }
679  else {
680  outmat.open((m_pathtxt+filename).c_str());
681  if(outmat.fail())
682  return StatusCode::FAILURE;
683  outmat.setf(std::ios::fixed);
684  outmat.setf(std::ios::showpoint);
685  outmat.precision(3);
686  outmat << "msizz: " << std::setw(6) << msizz << std::endl;
687  outmat << "nelem: " << std::setw(6) << nelem << std::endl;
688  outmat << "AlSpaMat version: " << std::setw(6) << version << std::endl;
689  }
690 
691  double melem=0;
692  long int ii, jj;
693  int32_t i, j;
694 
695  for (const datamap::value_type& p : m_ptr_map) {
696  melem = p.second;
697  elem(p.first, ii, jj); i=ii; j=jj; // just a type conversion
698  if(binary) {
699  outmat.write(charAddress((i)), sizeof (i));
700  outmat.write(charAddress((j)), sizeof (j));
701  outmat.write(charAddress((melem)), sizeof (melem));
702  }
703  else
704  outmat << std::setw(6) << i << std::setw(6) << j << std::setw(18) << melem << std::endl;
705  }
706  outmat.close();
707 
708  return StatusCode::SUCCESS;
709 }

Member Data Documentation

◆ m_matrix_type

int Trk::AlSymMatBase::m_matrix_type = 0
protectedinherited

Definition at line 100 of file AlSymMatBase.h.

◆ m_nele

long int Trk::AlSymMatBase::m_nele = 0
protectedinherited

Definition at line 104 of file AlSymMatBase.h.

◆ m_pathbin

std::string Trk::AlSpaMat::m_pathbin
protected

Definition at line 85 of file AlSpaMat.h.

◆ m_pathtxt

std::string Trk::AlSpaMat::m_pathtxt
protected

Definition at line 86 of file AlSpaMat.h.

◆ m_ptr_col

int* Trk::AlSpaMat::m_ptr_col = nullptr
protected

Definition at line 89 of file AlSpaMat.h.

◆ m_ptr_map

datamap Trk::AlSymMatBase::m_ptr_map
protectedinherited

Definition at line 101 of file AlSymMatBase.h.

◆ m_ptr_row

int* Trk::AlSpaMat::m_ptr_row = nullptr
protected

Definition at line 88 of file AlSpaMat.h.

◆ m_size

long int Trk::AlSymMatBase::m_size = 0
protectedinherited

Definition at line 103 of file AlSymMatBase.h.


The documentation for this class was generated from the following files:
AllowedVariables::e
e
Definition: AsgElectronSelectorTool.cxx:37
beamspotman.r
def r
Definition: beamspotman.py:674
Trk::AlSpaMat::m_ptr_row
int * m_ptr_row
Definition: AlSpaMat.h:88
Trk::const_mapiterator
datamap::const_iterator const_mapiterator
Definition: AlSymMatBase.h:30
athena.path
path
python interpreter configuration --------------------------------------—
Definition: athena.py:128
ClusterSeg::residual
@ residual
Definition: ClusterNtuple.h:20
index
Definition: index.py:1
hist_file_dump.d
d
Definition: hist_file_dump.py:142
Trk::AlSpaMat::elemr
virtual double & elemr(long int, long int) override final
Definition: AlSpaMat.cxx:143
Trk::indices
std::pair< long int, long int > indices
Definition: AlSymMatBase.h:24
Trk::AlSpaMat::copy
void copy(const AlSpaMat &m)
Definition: AlSpaMat.cxx:92
Trk::AlSymMatBase::m_ptr_map
datamap m_ptr_map
Definition: AlSymMatBase.h:101
Trk::AlSpaMat::elemc
virtual double elemc(long int, long int) const override final
Definition: AlSpaMat.cxx:179
keylayer_zslicemap.row
row
Definition: keylayer_zslicemap.py:155
UploadAMITag.l
list l
Definition: UploadAMITag.larcaf.py:157
JetTiledMap::N
@ N
Definition: TiledEtaPhiMap.h:44
Trk::AlSpaMat::m_pathbin
std::string m_pathbin
Definition: AlSpaMat.h:85
Trk::AlSpaMat::CheckMatVersion
virtual StatusCode CheckMatVersion(const std::string &, bool &) override final
Definition: AlSpaMat.cxx:712
atlasStyleMacro.icol
int icol
Definition: atlasStyleMacro.py:13
Trk::AlSpaMat::elem
void elem(const indices &, long int &, long int &) const
Definition: AlSpaMat.h:94
Trk::AlSymMatBase::AlSymMatBase
AlSymMatBase()
Definition: AlSymMatBase.cxx:20
TMatrixDSparse
class TMatrixTSparse< double > TMatrixDSparse
Definition: AlSpaMat.h:14
Trk::AlSymMatBase::size
long int size() const
Definition: AlSymMatBase.h:145
python.utils.AtlRunQueryDQUtils.p
p
Definition: AtlRunQueryDQUtils.py:209
lumiFormat.i
int i
Definition: lumiFormat.py:85
beamspotman.n
n
Definition: beamspotman.py:729
Trk::AlSpaMat::AlSpaMat
AlSpaMat()
Definition: AlSpaMat.cxx:43
Trk::AlSymMatBase::m_nele
long int m_nele
Definition: AlSymMatBase.h:104
Trk::AlSpaMat::RemoveDoF
void RemoveDoF(int, int nelem=1)
Definition: AlSpaMat.cxx:493
xAOD::double
double
Definition: CompositeParticle_v1.cxx:159
plotBeamSpotMon.b
b
Definition: plotBeamSpotMon.py:76
Trk::mapiterator
datamap::iterator mapiterator
Definition: AlSymMatBase.h:29
python.LumiBlobConversion.pos
pos
Definition: LumiBlobConversion.py:16
get_generator_info.version
version
Definition: get_generator_info.py:33
DeMoScan.index
string index
Definition: DeMoScan.py:362
a
TList * a
Definition: liststreamerinfos.cxx:10
Pythia8_RapidityOrderMPI.val
val
Definition: Pythia8_RapidityOrderMPI.py:14
CaloCellTimeCorrFiller.filename
filename
Definition: CaloCellTimeCorrFiller.py:23
Trk::AlSpaMat::m_ptr_col
int * m_ptr_col
Definition: AlSpaMat.h:89
Trk::AlSymMatBase::m_matrix_type
int m_matrix_type
Definition: AlSymMatBase.h:100
Trk::x
@ x
Definition: ParamDefs.h:55
test_pyathena.counter
counter
Definition: test_pyathena.py:15
Trk::AlSpaMat::m_pathtxt
std::string m_pathtxt
Definition: AlSpaMat.h:86
Trk::AlSymMatBase::m_size
long int m_size
Definition: AlSymMatBase.h:103
Trk::v
@ v
Definition: ParamDefs.h:78
fitman.k
k
Definition: fitman.py:528
python.SystemOfUnits.m
float m
Definition: SystemOfUnits.py:106
mapkey::key
key
Definition: TElectronEfficiencyCorrectionTool.cxx:37