MOCHA 0.9
|
00001 // Copyright 2009 by Jesus De Loera, David Haws, Jon Lee, Allison O'Hair, University 00002 // of California at Davis, IBM, Inc. All rights reserved. 00003 // 00004 // Distributed under the Eclipse Public License v 1.0. See ../LICENSE 00005 00006 // $Rev$ $Date$ 00007 00008 // matrix.h. Author David Haws 00009 00010 #ifndef MATRIX_H 00011 #define MATRIX_H 00012 00013 #ifdef HAVE_CONFIG_H 00014 #include "../config.h" 00015 #endif 00016 00017 #include <iostream> 00018 #include <set> 00019 00020 /* cannot have gmpxx without gmp */ 00021 #ifndef HAVE_LIBGMP 00022 #ifdef HAVE_LIBGMPXX 00023 #undef HAVE_LIBGMPXX 00024 #endif 00025 #endif 00026 00027 #ifdef HAVE_LIBGMPXX 00028 #include <gmpxx.h> 00029 #endif 00030 00031 using namespace::std; 00032 00033 class Matrix 00034 { 00035 public: 00036 Matrix(); // Does not allocate Entries. Sets rows = cols = 0 00037 Matrix(unsigned,unsigned); // Creates mxn matrix 00038 Matrix(unsigned); // Creates square matrix 00039 Matrix(unsigned,unsigned,int, int); // Creates mxn matrix filled with random integers ranged x to y 00040 ~Matrix(); 00041 Matrix(const Matrix& m); // Copy constructor 00042 Matrix& operator= (const Matrix& m); // Assignment operator 00043 00044 bool operator== (const Matrix &m); // Equal operator 00045 00046 // Indexed 0 - (rows - 1), 0 - (cols - 1) 00047 double& operator() (unsigned row, unsigned col); 00048 double operator() (unsigned row, unsigned col) const; 00049 Matrix operator* (const Matrix &B) const; 00050 Matrix operator* (const double &rhs) const; 00051 Matrix operator+ (const Matrix &B) const; 00052 Matrix operator- (const Matrix &B) const; 00053 friend std::ostream& operator<< (std::ostream& o, const Matrix &someMatrix); 00054 friend std::istream& operator>> (std::istream& in, Matrix &someMatrix); 00055 00056 // This returns the transpose of this matrix 00057 Matrix transpose(); 00058 00059 // This will output the matrix in matlab format with label 00060 void writeMatlab(std::ostream &o, string label); 00061 00062 //Matrix *subColumns(set<unsigned> &colSet); 00063 // Returns a submatrix using colSet as columns to include 00064 Matrix subColumns(const set<unsigned> &colSet); 00065 00066 //Matrix *subColumns(set<unsigned> &colSet); 00067 // Returns a submatrix using NOT colSet as columns to include 00068 Matrix subColumnsDiff(const set<unsigned> &colSet); 00069 00070 //given a mxn matrix, returns a mx1 matrix which is the sum of each row 00071 Matrix rowSum(); 00072 00073 // Returns the cofactor matrix determined by i and j. 00074 // rows and cols > 1 00075 Matrix cofactor(unsigned i, unsigned j); 00076 00077 // Returns the result of Gaussian Elimination 00078 //Matrix *GE (); 00079 Matrix GE (); 00080 00081 #ifdef HAVE_LIBGMPXX 00082 // Returns the rank of this matrix. Uses arbitrary precision and expects 00083 // all entries to be integer. Also will set rank of the matrix 00084 int GE_gmp (); 00085 int GE_gmp_verbose (); 00086 00087 #endif 00088 int rank (); // Returns the rank of the matrix 00089 00090 double trace(); // Computes the trace 00091 00092 long double det(); // Computes the determinant 00093 long double detCofactor(); // Computes the determinant using cofactor expansion 00094 long double detLU(); // Computes the determinant using SVD decomposition 00095 long double log_abs_det(); // Computes log(abs(det())); 00096 00097 // This uses dgesvd function of lapack to compute *this = U*S*V.transpose() 00098 void SVD(Matrix &U, Matrix &S, Matrix &V); 00099 00100 // This uses dgetrf function of lapack to compute *this = P*L*U. 00101 // returns the sign determined by the number of row transpositions in P 00102 double LU(Matrix &P, Matrix &L, Matrix &U); 00103 00104 // This uses the SVD decomposition function above to compute the inverse. 00105 Matrix inverse(); 00106 00107 // Uses LAPACK and calls dgesvd_. Counts number of singular values 00108 // that are not < 10^-16 * SIGMA_MAX 00109 int rank_LAPACK (); 00110 const int getCols (); 00111 const int getRows (); 00112 //double det(); 00113 // Swaps columns i and j 00114 void swapCols(unsigned i, unsigned j); 00115 // Swaps rows i and j 00116 void swapRows(unsigned i, unsigned j); 00117 static int printPadLength; 00118 00119 // If this matrix is a column vector, return the 2-norm. Else exit with error 00120 double twoNorm (); 00121 double twoNormSquared (); 00122 unsigned rows; 00123 unsigned cols; 00124 private: 00125 double **Entries; 00126 void deleteEntries (); 00127 int matrixRank; 00128 void initialize(unsigned,unsigned); 00129 }; 00130 00131 struct ltcolvec 00132 { 00133 bool operator()(const Matrix &m1, const Matrix &m2) const 00134 { 00135 int allEqual = 1; 00136 if ((m1.rows) < (m2.rows)) 00137 { 00138 return 1; 00139 } 00140 for(unsigned i=0;i<m1.rows;i++) 00141 { 00142 if ((m1)(i,0) != (m2)(i,0)) 00143 { 00144 allEqual = 0; 00145 } 00146 if ((m1)(i,0) < (m2)(i,0)) 00147 { 00148 return 1; 00149 } 00150 else if ((m1)(i,0) > (m2)(i,0)) 00151 { 00152 return 0; 00153 } 00154 } 00155 if (allEqual == 1) 00156 { 00157 return 0; 00158 } 00159 return 1; 00160 } 00161 }; 00162 00163 #endif 00164