jiri/LA_library
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

*** empty log message ***

Browse Source
This commit is contained in:
jiri 2005-10-12 11:14:55 +00:00
parent 25f8b1eb6c
commit 03bed24414
3 changed files with 7 additions and 8 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

9
mat.h
View File

@ -86,11 +86,6 @@ public:
const NRMat submatrix(const int fromrow, const int torow, const int fromcol, const int tocol) const; //there is also independent less efficient routine for generally indexed submatrix const NRMat submatrix(const int fromrow, const int torow, const int fromcol, const int tocol) const; //there is also independent less efficient routine for generally indexed submatrix
void gemm(const T &beta, const NRMat &a, const char transa, const NRMat &b, void gemm(const T &beta, const NRMat &a, const char transa, const NRMat &b,
const char transb, const T &alpha);//this = alpha*op( A )*op( B ) + beta*this const char transb, const T &alpha);//this = alpha*op( A )*op( B ) + beta*this
/*
void strassen(const T beta, const NRMat &a, const char transa, const NRMat &b,
const char transb, const T alpha);//this := alpha*op( A )*op( B ) + beta*this
void s_cutoff(const int,const int,const int,const int) const;
*/
void fprintf(FILE *f, const char *format, const int modulo) const; void fprintf(FILE *f, const char *format, const int modulo) const;
void fscanf(FILE *f, const char *format); void fscanf(FILE *f, const char *format);
const double norm(const T scalar=(T)0) const; const double norm(const T scalar=(T)0) const;
@ -105,11 +100,11 @@ public:
void gemm(const T &beta, const SparseMat<T> &a, const char transa, const NRMat &b, const char transb, const T &alpha);//this = alpha*op( A )*op( B ) + beta*this void gemm(const T &beta, const SparseMat<T> &a, const char transa, const NRMat &b, const char transb, const T &alpha);//this = alpha*op( A )*op( B ) + beta*this
inline void simplify() {}; //just for compatibility with sparse ones inline void simplify() {}; //just for compatibility with sparse ones
bool issymmetric() const {return 0;}; bool issymmetric() const {return 0;};
#ifndef NO_STRASSEN
//Strassen's multiplication (better than n^3, analogous syntax to gemm) //Strassen's multiplication (better than n^3, analogous syntax to gemm)
void strassen(const T beta, const NRMat &a, const char transa, const NRMat &b, const char transb, const T alpha);//this := alpha*op( A )*op( B ) + beta*this void strassen(const T beta, const NRMat &a, const char transa, const NRMat &b, const char transb, const T alpha);//this := alpha*op( A )*op( B ) + beta*this
void s_cutoff(const int,const int,const int,const int) const; void s_cutoff(const int,const int,const int,const int) const;
#endif
}; };
//due to mutual includes this has to be after full class declaration //due to mutual includes this has to be after full class declaration

3
strassen.cc
View File

@ -1,6 +1,7 @@
#include "mat.h" #include "mat.h"
#ifndef NO_STRASSEN
/*Strassen algorithm*/ /*Strassen algorithm*/
// called routine is fortran-compatible // called routine is fortran-compatible
extern "C" void fmm(const char c_transa,const char c_transb,const int m,const int n,const int k,const double alpha, extern "C" void fmm(const char c_transa,const char c_transb,const int m,const int n,const int k,const double alpha,
@ -23,4 +24,4 @@ copyonwrite();
//swap transpositions and order of matrices //swap transpositions and order of matrices
fmm(transb,transa,mm,nn,k,alpha,b,b.mm, a, a.mm, beta,*this, mm,NULL,0); fmm(transb,transa,mm,nn,k,alpha,b,b.mm, a, a.mm, beta,*this, mm,NULL,0);
} }
#endif

3
t.cc
View File

@ -194,12 +194,14 @@ cout <<" big1*big2 "<<big3[0][0]<<" time "<<clock()/((double) (CLOCKS_PER_SEC))-
*/ */
#ifndef NO_STRASSEN
if(0) if(0)
{ {
NRMat<double> atest, btest,ctest; NRMat<double> atest, btest,ctest;
{ {
int cc,c1,c2,c3; int cc,c1,c2,c3;
cin >>cc>>c1>>c2>>c3; cin >>cc>>c1>>c2>>c3;
atest.s_cutoff(cc,c1,c2,c3); atest.s_cutoff(cc,c1,c2,c3);
} }
cin>>atest; cin>>atest;
@ -237,6 +239,7 @@ for (int c=64; c<=512;c+=64)
cout <<"cutoff "<<c<<" big1*big2 "<<big4[0][0]<<" time "<<clock()/((double) (CLOCKS_PER_SEC))-t <<"\n"; cout <<"cutoff "<<c<<" big1*big2 "<<big4[0][0]<<" time "<<clock()/((double) (CLOCKS_PER_SEC))-t <<"\n";
} }
} }
#endif
if(0) if(0)
{ {