Skip to content
Snippets Groups Projects
Commit ca0393da authored by Jesper Andersson's avatar Jesper Andersson
Browse files

BFSparseMatrix now templated to allow double or float representation

parent 5fdaea2f
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
bfmatrix.cpp
+ 100
− 282
View file @ ca0393da
......@@ -17,216 +17,6 @@
namespace MISCMATHS {
//
// Member functions for SparseBFMatrix
//
//
// Transpose
//
boost::shared_ptr<BFMatrix> SparseBFMatrix::Transpose()
const
{
boost::shared_ptr<SparseBFMatrix> tm(new SparseBFMatrix(mp->t()));
return(tm);
}
//
// Concatenation of two matrices returning a third
//
void SparseBFMatrix::HorConcat(const BFMatrix& B, BFMatrix& AB) const
{
if (B.Nrows() && Nrows() != B.Nrows()) {throw BFMatrixException("SparseBFMatrix::HorConcat: Matrices must have same # of rows");}
try {
SparseBFMatrix& lAB = dynamic_cast<SparseBFMatrix&>(AB); // Means that output is sparse
lAB = *this;
lAB.HorConcat2MyRight(B);
}
catch (std::bad_cast) {
try {
FullBFMatrix& lAB = dynamic_cast<FullBFMatrix&>(AB); // Means that output is full
lAB = FullBFMatrix(this->AsMatrix());
lAB.HorConcat2MyRight(B);
}
catch (std::bad_cast) {
throw BFMatrixException("SparseBFMatrix::HorConcat: dynamic cast error");
}
}
}
void SparseBFMatrix::HorConcat(const NEWMAT::Matrix& B, BFMatrix& AB) const
{
if (B.Nrows() && int(Nrows()) != B.Nrows()) {throw BFMatrixException("SparseBFMatrix::HorConcat: Matrices must have same # of rows");}
try {
SparseBFMatrix& lAB = dynamic_cast<SparseBFMatrix&>(AB); // Means that output is sparse
lAB = *this;
lAB.HorConcat2MyRight(B);
}
catch (std::bad_cast) {
try {
FullBFMatrix& lAB = dynamic_cast<FullBFMatrix&>(AB); // Means that output is full
lAB = FullBFMatrix(this->AsMatrix());
lAB.HorConcat2MyRight(B);
}
catch (std::bad_cast) {
throw BFMatrixException("SparseBFMatrix::HorConcat: dynamic cast error");
}
}
}
void SparseBFMatrix::VertConcat(const BFMatrix& B, BFMatrix& AB) const
{
if (B.Ncols() && Ncols() != B.Ncols()) {throw BFMatrixException("SparseBFMatrix::VertConcat: Matrices must have same # of columns");}
try {
SparseBFMatrix& lAB = dynamic_cast<SparseBFMatrix&>(AB); // Means that output is sparse
lAB = *this;
lAB.VertConcatBelowMe(B);
}
catch (std::bad_cast) {
try {
FullBFMatrix& lAB = dynamic_cast<FullBFMatrix&>(AB); // Means that output is full
lAB = FullBFMatrix(this->AsMatrix());
lAB.VertConcatBelowMe(B);
}
catch (std::bad_cast) {
throw BFMatrixException("SparseBFMatrix::VertConcat: dynamic cast error");
}
}
}
void SparseBFMatrix::VertConcat(const NEWMAT::Matrix& B, BFMatrix& AB) const
{
if (B.Ncols() && int(Ncols()) != B.Ncols()) {throw BFMatrixException("SparseBFMatrix::VertConcat: Matrices must have same # of columns");}
try {
SparseBFMatrix& lAB = dynamic_cast<SparseBFMatrix&>(AB); // Means that output is sparse
lAB = *this;
lAB.VertConcatBelowMe(B);
}
catch (std::bad_cast) {
try {
FullBFMatrix& lAB = dynamic_cast<FullBFMatrix&>(AB); // Means that output is full
lAB = FullBFMatrix(this->AsMatrix());
lAB.VertConcatBelowMe(B);
}
catch (std::bad_cast) {
throw BFMatrixException("SparseBFMatrix::VertConcat: dynamic cast error");
}
}
}
//
// Concatenate another matrix to *this
//
void SparseBFMatrix::HorConcat2MyRight(const BFMatrix& B)
{
if (!B.Nrows()) return;
if (Nrows() != B.Nrows()) {throw BFMatrixException("SparseBFMatrix::HorConcat2MyRight: Matrices must have same # of rows");}
try {
const SparseBFMatrix& lB = dynamic_cast<const SparseBFMatrix&>(B); // Means that we want to concatenate a sparse matrix
*mp |= *(lB.mp);
}
catch (std::bad_cast) {
try {
const FullBFMatrix& lB = dynamic_cast<const FullBFMatrix&>(B); // Means that we want to concatenate a full
this->HorConcat2MyRight(lB.AsMatrix());
}
catch (std::bad_cast) {
throw BFMatrixException("SparseBFMatrix::HorConcat2MyRight: dynamic cast error");
}
}
}
void SparseBFMatrix::HorConcat2MyRight(const NEWMAT::Matrix& B)
{
if (!B.Nrows()) return;
if (int(Nrows()) != B.Nrows()) {throw BFMatrixException("SparseBFMatrix::HorConcat2MyRight: Matrices must have same # of rows");}
*mp |= B;
}
void SparseBFMatrix::VertConcatBelowMe(const BFMatrix& B)
{
if (!B.Ncols()) return;
if (Ncols() != B.Ncols()) {throw BFMatrixException("SparseBFMatrix::VertConcatBelowMe: Matrices must have same # of columns");}
try {
const SparseBFMatrix& lB = dynamic_cast<const SparseBFMatrix&>(B); // Means that we want to concatenate a sparse matrix
*mp &= *(lB.mp);
}
catch (std::bad_cast) {
try {
const FullBFMatrix& lB = dynamic_cast<const FullBFMatrix&>(B); // Means that we want to concatenate a full
this->VertConcatBelowMe(lB.AsMatrix());
}
catch (std::bad_cast) {
throw BFMatrixException("SparseBFMatrix::HorConcat2MyRight: dynamic cast error");
}
}
}
void SparseBFMatrix::VertConcatBelowMe(const NEWMAT::Matrix& B)
{
if (!B.Ncols()) return;
if (int(Ncols()) != B.Ncols()) {throw BFMatrixException("SparseBFMatrix::VertConcatBelowMe: Matrices must have same # of columns");}
*mp &= B;
}
// Multiply by vector
NEWMAT::ReturnMatrix SparseBFMatrix::MulByVec(const NEWMAT::ColumnVector& invec) const
{
if (invec.Nrows() != int(Ncols())) {throw BFMatrixException("Matrix-vector size mismatch");}
NEWMAT::ColumnVector outvec = *mp * invec;
outvec.Release();
return(outvec);
}
// Add another matrix to this one
void SparseBFMatrix::AddToMe(const BFMatrix& M, double s)
{
if (Ncols() != M.Ncols() || Nrows() != M.Nrows()) {
throw BFMatrixException("SparseBFMatrix::AddToMe: Matrix size mismatch");
}
try {
const SparseBFMatrix& lM = dynamic_cast<const SparseBFMatrix&>(M); // Add sparse matrix to this sparse matrix
if (s == 1.0) *mp += *(lM.mp);
else *mp += s * *(lM.mp);
}
catch (std::bad_cast) {
try {
const FullBFMatrix& lM = dynamic_cast<const FullBFMatrix&>(M); // Add full matrix to this sparse matrix
if (s == 1.0) *mp += SpMat<double>(lM.ReadAsMatrix());
else *mp += s * SpMat<double>(lM.ReadAsMatrix());
}
catch (std::bad_cast) {
throw BFMatrixException("SparseBFMatrix::AddToMe: dynamic cast error");
}
}
}
// Given A*x=b, solve for x
NEWMAT::ReturnMatrix SparseBFMatrix::SolveForx(const NEWMAT::ColumnVector& b,
MISCMATHS::MatrixType type,
double tol,
int miter) const
{
if (b.Nrows() != int(Nrows())) {
throw BFMatrixException("SparseBFMatrix::SolveForx: Matrix-vector size mismatch");
}
NEWMAT::ColumnVector x = mp->SolveForx(b,type,tol,miter);
x.Release();
return(x);
}
//
// Member functions for FullBFMatrix
//
......@@ -252,19 +42,24 @@ void FullBFMatrix::HorConcat(const BFMatrix& B, BFMatrix& AB) const
{
if (B.Nrows() && Nrows() != B.Nrows()) {throw BFMatrixException("FullBFMatrix::HorConcat: Matrices must have same # of rows");}
try {
FullBFMatrix& lAB = dynamic_cast<FullBFMatrix&>(AB); // This means output is a full matrix
lAB = *this;
lAB.HorConcat2MyRight(B);
FullBFMatrix *pAB = dynamic_cast<FullBFMatrix *>(&AB);
if (pAB) { // This means output is a full matrix
*pAB = *this;
pAB->HorConcat2MyRight(B);
}
catch (std::bad_cast) {
try {
SparseBFMatrix& lAB = dynamic_cast<SparseBFMatrix&>(AB);
lAB = SparseBFMatrix(this->AsMatrix());
lAB.HorConcat2MyRight(B);
else {
SparseBFMatrix<double> *psdAB = dynamic_cast<SparseBFMatrix<double> *>(&AB);
if (psdAB) {
*psdAB = SparseBFMatrix<double>(this->AsMatrix());
psdAB->HorConcat2MyRight(B);
}
catch (std::bad_cast) {
throw BFMatrixException("FullBFMatrix::HorConcat: dynamic cast error");
else {
SparseBFMatrix<float> *psfAB = dynamic_cast<SparseBFMatrix<float> *>(&AB);
if (psfAB) {
*psfAB = SparseBFMatrix<float>(this->AsMatrix());
psfAB->HorConcat2MyRight(B);
}
else throw BFMatrixException("FullBFMatrix::HorConcat: dynamic cast error");
}
}
}
......@@ -273,19 +68,24 @@ void FullBFMatrix::HorConcat(const NEWMAT::Matrix& B, BFMatrix& AB) const
{
if (B.Nrows() && int(Nrows()) != B.Nrows()) {throw BFMatrixException("FullBFMatrix::HorConcat: Matrices must have same # of rows");}
try {
FullBFMatrix& lAB = dynamic_cast<FullBFMatrix&>(AB); // This means output is a full matrix
lAB = *this;
lAB.HorConcat2MyRight(B);
FullBFMatrix *pAB = dynamic_cast<FullBFMatrix *>(&AB);
if (pAB) { // This means output is a full matrix
*pAB = *this;
pAB->HorConcat2MyRight(B);
}
catch (std::bad_cast) {
try {
SparseBFMatrix& lAB = dynamic_cast<SparseBFMatrix&>(AB);
lAB = SparseBFMatrix(this->AsMatrix());
lAB.HorConcat2MyRight(B);
else {
SparseBFMatrix<double> *psdAB = dynamic_cast<SparseBFMatrix<double> *>(&AB);
if (psdAB) {
*psdAB = SparseBFMatrix<double>(this->AsMatrix());
psdAB->HorConcat2MyRight(B);
}
catch (std::bad_cast) {
throw BFMatrixException("FullBFMatrix::HorConcat: dynamic cast error");
else {
SparseBFMatrix<float> *psfAB = dynamic_cast<SparseBFMatrix<float> *>(&AB);
if (psfAB) {
*psfAB = SparseBFMatrix<float>(this->AsMatrix());
psfAB->HorConcat2MyRight(B);
}
else throw BFMatrixException("FullBFMatrix::HorConcat: dynamic cast error");
}
}
}
......@@ -294,19 +94,24 @@ void FullBFMatrix::VertConcat(const BFMatrix& B, BFMatrix& AB) const
{
if (B.Ncols() && Ncols() != B.Ncols()) {throw BFMatrixException("FullBFMatrix::VertConcat: Matrices must have same # of columns");}
try {
FullBFMatrix& lAB = dynamic_cast<FullBFMatrix&>(AB); // This means output is a full matrix
lAB = *this;
lAB.VertConcatBelowMe(B);
FullBFMatrix *pAB = dynamic_cast<FullBFMatrix *>(&AB);
if (pAB) { // This means output is a full matrix
*pAB = *this;
pAB->VertConcatBelowMe(B);
}
catch (std::bad_cast) {
try {
SparseBFMatrix& lAB = dynamic_cast<SparseBFMatrix&>(AB);
lAB = SparseBFMatrix(this->AsMatrix());
lAB.VertConcatBelowMe(B);
else {
SparseBFMatrix<double> *psdAB = dynamic_cast<SparseBFMatrix<double> *>(&AB);
if (psdAB) {
*psdAB = SparseBFMatrix<double>(this->AsMatrix());
psdAB->VertConcatBelowMe(B);
}
catch (std::bad_cast) {
throw BFMatrixException("FullBFMatrix::VertConcat: dynamic cast error");
else {
SparseBFMatrix<float> *psfAB = dynamic_cast<SparseBFMatrix<float> *>(&AB);
if (psfAB) {
*psfAB = SparseBFMatrix<float>(this->AsMatrix());
psfAB->VertConcatBelowMe(B);
}
else throw BFMatrixException("FullBFMatrix::VertConcat: dynamic cast error");
}
}
}
......@@ -315,19 +120,24 @@ void FullBFMatrix::VertConcat(const NEWMAT::Matrix& B, BFMatrix& AB) const
{
if (B.Ncols() && int(Ncols()) != B.Ncols()) {throw BFMatrixException("FullBFMatrix::VertConcat: Matrices must have same # of columns");}
try {
FullBFMatrix& lAB = dynamic_cast<FullBFMatrix&>(AB); // This means output is a full matrix
lAB = *this;
lAB.VertConcatBelowMe(B);
FullBFMatrix *pAB = dynamic_cast<FullBFMatrix *>(&AB);
if (pAB) { // This means output is a full matrix
*pAB = *this;
pAB->VertConcatBelowMe(B);
}
catch (std::bad_cast) {
try {
SparseBFMatrix& lAB = dynamic_cast<SparseBFMatrix&>(AB);
lAB = SparseBFMatrix(this->AsMatrix());
lAB.VertConcatBelowMe(B);
else {
SparseBFMatrix<double> *psdAB = dynamic_cast<SparseBFMatrix<double> *>(&AB);
if (psdAB) {
*psdAB = SparseBFMatrix<double>(this->AsMatrix());
psdAB->VertConcatBelowMe(B);
}
catch (std::bad_cast) {
throw BFMatrixException("FullBFMatrix::VertConcat: dynamic cast error");
else {
SparseBFMatrix<float> *psfAB = dynamic_cast<SparseBFMatrix<float> *>(&AB);
if (psfAB) {
*psfAB = SparseBFMatrix<float>(this->AsMatrix());
psfAB->VertConcatBelowMe(B);
}
else throw BFMatrixException("FullBFMatrix::VertConcat: dynamic cast error");
}
}
}
......@@ -341,17 +151,21 @@ void FullBFMatrix::HorConcat2MyRight(const BFMatrix& B)
if (Nrows() != B.Nrows()) {throw BFMatrixException("FullBFMatrix::HorConcat2MyRight: Matrices must have same # of rows");}
try {
const FullBFMatrix& lB = dynamic_cast<const FullBFMatrix&>(B);
*mp |= *(lB.mp);
const FullBFMatrix *pB = dynamic_cast<const FullBFMatrix *>(&B);
if (pB) { // If B was full
*mp |= *(pB->mp);
}
catch (std::bad_cast) {
try {
const SparseBFMatrix& lB = dynamic_cast<const SparseBFMatrix&>(B);
this->HorConcat2MyRight(lB.AsMatrix());
else {
const SparseBFMatrix<double> *psdB = dynamic_cast<const SparseBFMatrix<double> *>(&B);
if (psdB) {
this->HorConcat2MyRight(psdB->AsMatrix());
}
catch (std::bad_cast) {
throw BFMatrixException("FullBFMatrix::HorConcat2MyRight: dynamic cast error");
else {
const SparseBFMatrix<float> *psfB = dynamic_cast<const SparseBFMatrix<float> *>(&B);
if (psfB) {
this->HorConcat2MyRight(psfB->AsMatrix());
}
else throw BFMatrixException("FullBFMatrix::HorConcat2MyRight: dynamic cast error");
}
}
}
......@@ -370,17 +184,21 @@ void FullBFMatrix::VertConcatBelowMe(const BFMatrix& B)
if (Ncols() != B.Ncols()) {throw BFMatrixException("FullBFMatrix::VertConcatBelowMe: Matrices must have same # of columns");}
try {
const FullBFMatrix& lB = dynamic_cast<const FullBFMatrix&>(B);
*mp &= *(lB.mp);
const FullBFMatrix *pB = dynamic_cast<const FullBFMatrix *>(&B);
if (pB) { // Means B is full
*mp &= *(pB->mp);
}
catch (std::bad_cast) {
try {
const SparseBFMatrix& lB = dynamic_cast<const SparseBFMatrix&>(B);
this->VertConcatBelowMe(lB.AsMatrix());
else {
const SparseBFMatrix<double> *psdB = dynamic_cast<const SparseBFMatrix<double> *>(&B);
if (psdB) {
this->VertConcatBelowMe(psdB->AsMatrix());
}
catch (std::bad_cast) {
throw BFMatrixException("FullBFMatrix::HorConcatBelowMe: dynamic cast error");
else {
const SparseBFMatrix<float> *psfB = dynamic_cast<const SparseBFMatrix<float> *>(&B);
if (psfB) {
this->VertConcatBelowMe(psfB->AsMatrix());
}
else throw BFMatrixException("FullBFMatrix::HorConcatBelowMe: dynamic cast error");
}
}
}
......@@ -417,17 +235,17 @@ void FullBFMatrix::AddToMe(const BFMatrix& m, double s)
throw BFMatrixException("FullBFMatrix::AddToMe: Matrix size mismatch");
}
try {
const FullBFMatrix& lm = dynamic_cast<const FullBFMatrix&>(m);
*mp += s*(*lm.mp);
const FullBFMatrix *pm = dynamic_cast<const FullBFMatrix *>(&m);
if (pm) { // If m is full matrix
*mp += s*(*(pm->mp));
}
catch (std::bad_cast) {
try {
const SparseBFMatrix& lm = dynamic_cast<const SparseBFMatrix&>(m);
*mp += s*lm.AsMatrix();
}
catch (std::bad_cast) {
throw BFMatrixException("FullBFMatrix::AddToMe: dynamic cast error");
else {
const SparseBFMatrix<double> *psdm = dynamic_cast<const SparseBFMatrix<double> *>(&m);
if (psdm) *mp += s*psdm->AsMatrix();
else {
const SparseBFMatrix<float> *psfm = dynamic_cast<const SparseBFMatrix<float> *>(&m);
if (psfm) *mp += s*psfm->AsMatrix();
else throw BFMatrixException("FullBFMatrix::AddToMe: dynamic cast error");
}
}
}
......
bfmatrix.h
+ 234
− 16
View file @ ca0393da
......@@ -49,6 +49,8 @@ public:
~BFMatrixException() throw() {}
};
enum BFMatrixPrecisionType {BFMatrixDoublePrecision, BFMatrixFloatPrecision};
class BFMatrix
{
protected:
......@@ -70,7 +72,8 @@ public:
virtual void Print(const std::string fname=std::string("")) const = 0;
// Setting, deleting or resizing the whole sparse matrix.
virtual void SetMatrix(const MISCMATHS::SpMat<double>& M) = 0;
// virtual void SetMatrix(const MISCMATHS::SpMat<double>& M) = 0;
// virtual void SetMatrix(const MISCMATHS::SpMat<float>& M) = 0;
virtual void SetMatrix(const NEWMAT::Matrix& M) = 0;
virtual void Clear() = 0;
virtual void Resize(unsigned int m, unsigned int n) = 0;
......@@ -114,26 +117,27 @@ public:
int miter=200) const = 0;
};
template<class T>
class SparseBFMatrix : public BFMatrix
{
private:
boost::shared_ptr<MISCMATHS::SpMat<double> > mp;
boost::shared_ptr<MISCMATHS::SpMat<T> > mp;
public:
// Constructors, destructor and assignment
SparseBFMatrix()
: mp(boost::shared_ptr<MISCMATHS::SpMat<double> >(new MISCMATHS::SpMat<double>())) {}
: mp(boost::shared_ptr<MISCMATHS::SpMat<T> >(new MISCMATHS::SpMat<T>())) {}
SparseBFMatrix(unsigned int m, unsigned int n)
: mp(boost::shared_ptr<MISCMATHS::SpMat<double> >(new MISCMATHS::SpMat<double>(m,n))) {}
: mp(boost::shared_ptr<MISCMATHS::SpMat<T> >(new MISCMATHS::SpMat<T>(m,n))) {}
SparseBFMatrix(unsigned int m, unsigned int n, const unsigned int *irp, const unsigned int *jcp, const double *sp)
: mp(boost::shared_ptr<MISCMATHS::SpMat<double> >(new MISCMATHS::SpMat<double>(m,n,irp,jcp,sp))) {}
SparseBFMatrix(const MISCMATHS::SpMat<double>& M)
: mp(boost::shared_ptr<MISCMATHS::SpMat<double> >(new MISCMATHS::SpMat<double>(M))) {}
: mp(boost::shared_ptr<MISCMATHS::SpMat<T> >(new MISCMATHS::SpMat<T>(m,n,irp,jcp,sp))) {}
SparseBFMatrix(const MISCMATHS::SpMat<T>& M)
: mp(boost::shared_ptr<MISCMATHS::SpMat<T> >(new MISCMATHS::SpMat<T>(M))) {}
SparseBFMatrix(const NEWMAT::Matrix& M)
: mp(boost::shared_ptr<MISCMATHS::SpMat<double> >(new MISCMATHS::SpMat<double>(M))) {}
: mp(boost::shared_ptr<MISCMATHS::SpMat<T> >(new MISCMATHS::SpMat<T>(M))) {}
virtual ~SparseBFMatrix() {}
virtual const SparseBFMatrix& operator=(const SparseBFMatrix& M) {
mp = boost::shared_ptr<MISCMATHS::SpMat<double> >(new MISCMATHS::SpMat<double>(*(M.mp))); return(*this);
virtual const SparseBFMatrix& operator=(const SparseBFMatrix<T>& M) {
mp = boost::shared_ptr<MISCMATHS::SpMat<T> >(new MISCMATHS::SpMat<T>(*(M.mp))); return(*this);
}
// Access as NEWMAT::Matrix
......@@ -147,11 +151,12 @@ public:
virtual void Print(const std::string fname=std::string("")) const {mp->Print(fname);}
// Setting, deleting or resizing the whole sparse matrix.
virtual void SetMatrix(const MISCMATHS::SpMat<double>& M) {mp = boost::shared_ptr<MISCMATHS::SpMat<double> >(new MISCMATHS::SpMat<double>(M));}
virtual void SetMatrix(const NEWMAT::Matrix& M) {mp = boost::shared_ptr<MISCMATHS::SpMat<double> >(new MISCMATHS::SpMat<double>(M));}
virtual void SetMatrixPtr(boost::shared_ptr<MISCMATHS::SpMat<double> >& mptr) {mp = mptr;}
virtual void Clear() {mp = boost::shared_ptr<MISCMATHS::SpMat<double> >(new MISCMATHS::SpMat<double>());}
virtual void Resize(unsigned int m, unsigned int n) {mp = boost::shared_ptr<MISCMATHS::SpMat<double> >(new MISCMATHS::SpMat<double>(m,n));}
virtual void SetMatrix(const MISCMATHS::SpMat<T>& M) {mp = boost::shared_ptr<MISCMATHS::SpMat<T> >(new MISCMATHS::SpMat<T>(M));}
// virtual void SetMatrix(const MISCMATHS::SpMat<float>& M) {mp = boost::shared_ptr<MISCMATHS::SpMat<float> >(new MISCMATHS::SpMat<float>(M));}
virtual void SetMatrix(const NEWMAT::Matrix& M) {mp = boost::shared_ptr<MISCMATHS::SpMat<T> >(new MISCMATHS::SpMat<T>(M));}
virtual void SetMatrixPtr(boost::shared_ptr<MISCMATHS::SpMat<T> >& mptr) {mp = mptr;}
virtual void Clear() {mp = boost::shared_ptr<MISCMATHS::SpMat<T> >(new MISCMATHS::SpMat<T>());}
virtual void Resize(unsigned int m, unsigned int n) {mp = boost::shared_ptr<MISCMATHS::SpMat<T> >(new MISCMATHS::SpMat<T>(m,n));}
// Accessing values
virtual double Peek(unsigned int r, unsigned int c) const {return(mp->Peek(r,c));}
......@@ -221,6 +226,7 @@ public:
// Setting, deleting or resizing the whole matrix.
virtual void SetMatrix(const MISCMATHS::SpMat<double>& M) {mp = boost::shared_ptr<NEWMAT::Matrix>(new NEWMAT::Matrix(M.AsNEWMAT()));}
virtual void SetMatrix(const MISCMATHS::SpMat<float>& M) {mp = boost::shared_ptr<NEWMAT::Matrix>(new NEWMAT::Matrix(M.AsNEWMAT()));}
virtual void SetMatrix(const NEWMAT::Matrix& M) {mp = boost::shared_ptr<NEWMAT::Matrix>(new NEWMAT::Matrix(M));}
virtual void SetMatrixPtr(boost::shared_ptr<NEWMAT::Matrix>& mptr) {mp = mptr;}
virtual void Clear() {mp->ReSize(0,0);}
......@@ -266,6 +272,218 @@ public:
};
//
// Here comes member functions for SparseBFMatrix. Since it is templated
// these need to go here rather than in bfmatrix.cpp.
//
//
// Member functions for SparseBFMatrix
//
//
// Transpose
//
template<class T>
boost::shared_ptr<BFMatrix> SparseBFMatrix<T>::Transpose() const
{
boost::shared_ptr<SparseBFMatrix<T> > tm(new SparseBFMatrix<T>(mp->t()));
return(tm);
}
//
// Concatenation of two matrices returning a third
//
template<class T>
void SparseBFMatrix<T>::HorConcat(const BFMatrix& B, BFMatrix& AB) const
{
if (B.Nrows() && Nrows() != B.Nrows()) {throw BFMatrixException("SparseBFMatrix::HorConcat: Matrices must have same # of rows");}
SparseBFMatrix<T> *pAB = dynamic_cast<SparseBFMatrix<T> *>(&AB);
if (pAB) { // Means that output is sparse of type T
*pAB = *this;
pAB->HorConcat2MyRight(B);
}
else {
FullBFMatrix *fpAB = dynamic_cast<FullBFMatrix *>(&AB);
if (fpAB) { // Means that output is full
*fpAB = FullBFMatrix(this->AsMatrix());
fpAB->HorConcat2MyRight(B);
}
else throw BFMatrixException("SparseBFMatrix::HorConcat: dynamic cast error");
}
}
template<class T>
void SparseBFMatrix<T>::HorConcat(const NEWMAT::Matrix& B, BFMatrix& AB) const
{
if (B.Nrows() && int(Nrows()) != B.Nrows()) {throw BFMatrixException("SparseBFMatrix::HorConcat: Matrices must have same # of rows");}
SparseBFMatrix<T> *pAB = dynamic_cast<SparseBFMatrix<T> *>(&AB);
if (pAB) { // Means that output is sparse
*pAB = *this;
pAB->HorConcat2MyRight(B);
}
else {
FullBFMatrix *fpAB = dynamic_cast<FullBFMatrix *>(&AB);
if (fpAB) {// Means that output is full
*fpAB = FullBFMatrix(this->AsMatrix());
fpAB->HorConcat2MyRight(B);
}
else throw BFMatrixException("SparseBFMatrix::HorConcat: dynamic cast error");
}
}
template<class T>
void SparseBFMatrix<T>::VertConcat(const BFMatrix& B, BFMatrix& AB) const
{
if (B.Ncols() && Ncols() != B.Ncols()) {throw BFMatrixException("SparseBFMatrix::VertConcat: Matrices must have same # of columns");}
SparseBFMatrix<T> *pAB = dynamic_cast<SparseBFMatrix<T> *>(&AB);
if (pAB) { // Means that output is sparse
*pAB = *this;
pAB->VertConcatBelowMe(B);
}
else {
FullBFMatrix *fpAB = dynamic_cast<FullBFMatrix *>(&AB);
if (fpAB) { // Means that output is full
*fpAB = FullBFMatrix(this->AsMatrix());
fpAB->VertConcatBelowMe(B);
}
else throw BFMatrixException("SparseBFMatrix::VertConcat: dynamic cast error");
}
}
template<class T>
void SparseBFMatrix<T>::VertConcat(const NEWMAT::Matrix& B, BFMatrix& AB) const
{
if (B.Ncols() && int(Ncols()) != B.Ncols()) {throw BFMatrixException("SparseBFMatrix::VertConcat: Matrices must have same # of columns");}
SparseBFMatrix<T> *pAB = dynamic_cast<SparseBFMatrix<T> *>(&AB);
if (pAB) { // Means that output is sparse
*pAB = *this;
pAB->VertConcatBelowMe(B);
}
else {
FullBFMatrix *fpAB = dynamic_cast<FullBFMatrix *>(&AB);
if (fpAB) { // Means that output is full
*fpAB = FullBFMatrix(this->AsMatrix());
fpAB->VertConcatBelowMe(B);
}
else throw BFMatrixException("SparseBFMatrix::VertConcat: dynamic cast error");
}
}
//
// Concatenate another matrix to *this
//
template<class T>
void SparseBFMatrix<T>::HorConcat2MyRight(const BFMatrix& B)
{
if (!B.Nrows()) return;
if (Nrows() != B.Nrows()) {throw BFMatrixException("SparseBFMatrix::HorConcat2MyRight: Matrices must have same # of rows");}
const SparseBFMatrix<T> *pB = dynamic_cast<const SparseBFMatrix<T> *>(&B);
if (pB) { // Means that we want to concatenate a sparse matrix
*mp |= *(pB->mp);
}
else {
const FullBFMatrix *fpB = dynamic_cast<const FullBFMatrix *>(&B);
if (fpB) { // Means that we want to concatenate a full
this->HorConcat2MyRight(fpB->AsMatrix());
}
else throw BFMatrixException("SparseBFMatrix::HorConcat2MyRight: dynamic cast error");
}
}
template<class T>
void SparseBFMatrix<T>::HorConcat2MyRight(const NEWMAT::Matrix& B)
{
if (!B.Nrows()) return;
if (int(Nrows()) != B.Nrows()) {throw BFMatrixException("SparseBFMatrix::HorConcat2MyRight: Matrices must have same # of rows");}
*mp |= B;
}
template<class T>
void SparseBFMatrix<T>::VertConcatBelowMe(const BFMatrix& B)
{
if (!B.Ncols()) return;
if (Ncols() != B.Ncols()) {throw BFMatrixException("SparseBFMatrix::VertConcatBelowMe: Matrices must have same # of columns");}
const SparseBFMatrix<T> *pB = dynamic_cast<const SparseBFMatrix<T> *>(&B);
if (pB) { // Means that we want to concatenate a sparse matrix
*mp &= *(pB->mp);
}
else {
const FullBFMatrix *fpB = dynamic_cast<const FullBFMatrix *>(&B);
if (fpB) { // Means that we want to concatenate a full
this->VertConcatBelowMe(fpB->AsMatrix());
}
else throw BFMatrixException("SparseBFMatrix::VertConcatBelowMe: dynamic cast error");
}
}
template<class T>
void SparseBFMatrix<T>::VertConcatBelowMe(const NEWMAT::Matrix& B)
{
if (!B.Ncols()) return;
if (int(Ncols()) != B.Ncols()) {throw BFMatrixException("SparseBFMatrix::VertConcatBelowMe: Matrices must have same # of columns");}
*mp &= B;
}
// Multiply by vector
template<class T>
NEWMAT::ReturnMatrix SparseBFMatrix<T>::MulByVec(const NEWMAT::ColumnVector& invec) const
{
if (invec.Nrows() != int(Ncols())) {throw BFMatrixException("Matrix-vector size mismatch");}
NEWMAT::ColumnVector outvec = *mp * invec;
outvec.Release();
return(outvec);
}
// Add another matrix to this one
template<class T>
void SparseBFMatrix<T>::AddToMe(const BFMatrix& M, double s)
{
if (Ncols() != M.Ncols() || Nrows() != M.Nrows()) {
throw BFMatrixException("SparseBFMatrix::AddToMe: Matrix size mismatch");
}
const SparseBFMatrix<T> *pM = dynamic_cast<const SparseBFMatrix<T> *>(&M);
if (pM) { // Add sparse matrix to this sparse matrix
if (s == 1.0) *mp += *(pM->mp);
else *mp += s * *(pM->mp);
}
else {
const FullBFMatrix *fpM = dynamic_cast<const FullBFMatrix *>(&M);
if (fpM) { // Add full matrix to this sparse matrix
if (s == 1.0) *mp += SpMat<T>(fpM->ReadAsMatrix());
else *mp += s * SpMat<T>(fpM->ReadAsMatrix());
}
else throw BFMatrixException("SparseBFMatrix::AddToMe: dynamic cast error");
}
}
// Given A*x=b, solve for x
template<class T>
NEWMAT::ReturnMatrix SparseBFMatrix<T>::SolveForx(const NEWMAT::ColumnVector& b,
MISCMATHS::MatrixType type,
double tol,
int miter) const
{
if (b.Nrows() != int(Nrows())) {
throw BFMatrixException("SparseBFMatrix::SolveForx: Matrix-vector size mismatch");
}
NEWMAT::ColumnVector x = mp->SolveForx(b,type,tol,miter);
x.Release();
return(x);
}
} // End namespace MISCMATHS
#endif // End #ifndef BFMatrix_h
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment