MatrixSparse.h

Go to the documentation of this file.

// Copyright 2019-2020 CERN and copyright holders of ALICE O2.
// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders.
// All rights not expressly granted are reserved.
//
// This software is distributed under the terms of the GNU General Public
// License v3 (GPL Version 3), copied verbatim in the file "COPYING".
//
// In applying this license CERN does not waive the privileges and immunities
// granted to it by virtue of its status as an Intergovernmental Organization
// or submit itself to any jurisdiction.
 
 
#ifndef ALICEO2_FWDALIGN_MATRIXSPARSE_H
#define ALICEO2_FWDALIGN_MATRIXSPARSE_H
 
#include "ForwardAlign/MatrixSq.h"
#include "ForwardAlign/VectorSparse.h"
 
namespace o2
{
namespace fwdalign
{
 
class MatrixSparse : public MatrixSq
{
 public:
  MatrixSparse() = default;
 
  MatrixSparse(Int_t size);
 
  MatrixSparse(const MatrixSparse& mat);
 
  ~MatrixSparse() override { Clear(); }
 
  VectorSparse* GetRow(Int_t ir) const { return (ir < fNcols) ? fVecs[ir] : nullptr; }
  VectorSparse* GetRowAdd(Int_t ir);
 
  Int_t GetSize() const override { return fNrows; }
  virtual Int_t GetNRows() const { return fNrows; }
  virtual Int_t GetNCols() const { return fNcols; }
 
  void Clear(Option_t* option = "") override;
  void Reset() override
  {
    for (int i = fNcols; i--;) {
      GetRow(i)->Reset();
    }
  }
  void Print(Option_t* option = "") const override;
  MatrixSparse& operator=(const MatrixSparse& src);
  Double_t& operator()(Int_t row, Int_t col) override;
  Double_t operator()(Int_t row, Int_t col) const override;
  void SetToZero(Int_t row, Int_t col);
 
  Float_t GetDensity() const override;
 
  Double_t DiagElem(Int_t r) const override;
  Double_t& DiagElem(Int_t r) override;
 
  void SortIndices(Bool_t valuesToo = kFALSE);
 
  void MultiplyByVec(const TVectorD& vecIn, TVectorD& vecOut) const override;
 
  void MultiplyByVec(const Double_t* vecIn, Double_t* vecOut) const override;
 
  void AddToRow(Int_t r, Double_t* valc, Int_t* indc, Int_t n) override;
 
 protected:
  VectorSparse** fVecs = nullptr;
 
  ClassDefOverride(MatrixSparse, 0);
};
 
//___________________________________________________
inline void MatrixSparse::MultiplyByVec(const TVectorD& vecIn, TVectorD& vecOut) const
{
  MultiplyByVec((Double_t*)vecIn.GetMatrixArray(), (Double_t*)vecOut.GetMatrixArray());
}
 
//___________________________________________________
inline void MatrixSparse::SetToZero(Int_t row, Int_t col)
{
  if (IsSymmetric() && col > row) {
    Swap(row, col);
  }
  VectorSparse* rowv = GetRow(row);
  if (rowv) {
    rowv->SetToZero(col);
  }
}
 
//___________________________________________________
inline Double_t MatrixSparse::operator()(Int_t row, Int_t col) const
{
  if (IsSymmetric() && col > row) {
    Swap(row, col);
  }
  VectorSparse* rowv = GetRow(row);
  if (!rowv) {
    return 0;
  }
  return rowv->FindIndex(col);
}
 
//___________________________________________________
inline Double_t& MatrixSparse::operator()(Int_t row, Int_t col)
{
  //  printf("M: findindexAdd\n");
  if (IsSymmetric() && col > row) {
    Swap(row, col);
  }
  VectorSparse* rowv = GetRowAdd(row);
  if (col >= fNcols) {
    fNcols = col + 1;
  }
  return rowv->FindIndexAdd(col);
}
 
//___________________________________________________
inline Double_t MatrixSparse::DiagElem(Int_t row) const
{
  VectorSparse* rowv = GetRow(row);
  if (!rowv) {
    return 0;
  }
  if (IsSymmetric()) {
    return (rowv->GetNElems() > 0 && rowv->GetLastIndex() == row) ? rowv->GetLastElem() : 0.;
  } else {
    return rowv->FindIndex(row);
  }
}
 
//___________________________________________________
inline Double_t& MatrixSparse::DiagElem(Int_t row)
{
  VectorSparse* rowv = GetRowAdd(row);
  if (row >= fNcols) {
    fNcols = row + 1;
  }
  if (IsSymmetric()) {
    return (rowv->GetNElems() > 0 && rowv->GetLastIndex() == row) ? rowv->GetLastElem() : rowv->FindIndexAdd(row);
  } else {
    return rowv->FindIndexAdd(row);
  }
}
 
} // namespace fwdalign
} // namespace o2
 
#endif