home *** CD-ROM | disk | FTP | other *** search
/ Power Programming / powerprogramming1994.iso / progtool / microcrn / issue_40.arc / DAIMS.ARC / MATRIX.HXX < prev    next >
Text File  |  1988-02-10  |  8KB  |  231 lines
  1. #include <stream.h>
  2. #define TRUE 1
  3. #define FALSE 0
  4. #define EOL 10            /*  end of line */
  5. #define TINY ((double)1e-20)
  6. #define NL "\n"
  7.  
  8. /* 
  9. -*++ class _Vector: thing to build matrices out of
  10. ** 
  11. ** (*++ history: 
  12. **      4 Dec 87    Bruce Eckel    Creation date
  13. ** ++*)
  14. ** 
  15. ** (*++ detailed: 
  16. ** ++*)
  17. */
  18.  
  19. class _Vector    {
  20.   double *v_d;
  21.   int    v_sz;
  22. public:
  23.   _Vector(int sz) { v_d = new double[v_sz = sz]; }
  24.   ~_Vector() { delete v_d; }
  25.   void error(int x) { cout << "vector index " << x << " out of bounds\n"; }
  26.   double& operator[](int x)  { 
  27.     if (x >= 0 && x < v_sz)
  28.       return v_d[x]; 
  29.     else { error(x); return v_d[x];}
  30.   }
  31. };
  32.  
  33. class Cheb_vector;
  34.  
  35. /* 
  36. -*++ class matrix: double matrix class.  Rectangular, any size.
  37. ** 
  38. ** (*++ history: 
  39. **      4 Dec 87    Bruce Eckel    Creation date.  Was actually created
  40. ** long before this, but the GNUemacs documentation system was just made.
  41. ** ++*)
  42. ** 
  43. ** (*++ detailed: This class can be used by itself or can be derived from.  It
  44. ** includes determinants, inversion, and several statistical functions as well
  45. ** as the usual matrix stuff.  You can index into a matrix using the m[][]
  46. ** construct; this does bounds checking for you.
  47. ** ++*)
  48. */
  49.  
  50. class matrix    {
  51.   struct matrep {
  52.     _Vector** m_vec;
  53.     int    m_rows;
  54.     int m_cols;
  55.     int n;            /* reference count */
  56.   };
  57.   matrep *p;
  58. public:
  59.   matrix(int rows = 1, int columns = 1, double initval = 0);
  60.   matrix(matrix& x);
  61.   matrix(char * initfile); /* initialize from a "standard" matrix file */
  62.   matrix(char * flag, int dimension); /* create an identity or random matrix */
  63.   ~matrix();
  64.   int rows() { return p->m_rows; };    /* max rows in matrix */
  65.   int cols() { return p->m_cols; };    /* max cols in matrix */
  66.   _Vector& operator[](int x) { 
  67.     if (x >= 0 && x < rows())
  68.       return *(p->m_vec[x]); 
  69.     else { 
  70.       cerr << "row index out of bounds: " << x << "\n"; 
  71.       return *(p->m_vec[x]);
  72.     }
  73.   }
  74.   void data() { cout << "rows = " << rows() << " cols = " << cols() 
  75.           << " ref count = " << p->n << "\n"; };
  76.   friend ostream& operator<<(ostream &s, matrix& m); 
  77.   friend istream& operator>>(istream &s, matrix& m); 
  78.   void write_standard(char * filename, char * msg = ""); 
  79.                   /* write matrix to `standard' file */
  80.   void file_get(char * infile);    /* read a matrix in from a file */
  81.   void file_put(char * outfile); /* write a matrix out to a file */
  82.   matrix & operator+(matrix&); /* add two matrices */
  83.   matrix & operator+(double arg); /* add a scalar to each el */
  84.   matrix & operator-(matrix& arg); /* subtract two matrices */
  85.   matrix & operator-(); /* unary minus */
  86.   matrix & operator*(double arg); /* multiply by a scalar */
  87.   matrix & operator*(matrix& arg); /* matrix multiply */
  88.   matrix & operator=(matrix&);    /* matrix assignment */
  89.   int operator==(matrix& arg); /* test for equivalence */
  90.   int operator!=(matrix& arg); /* test for inequivalence */
  91.   void set_row(int row, double value); /* set a row to a value */
  92.   void set_column(int column, double value); /* set a column to a value */
  93.   void switch_columns(int col1, int col2); /* swap two columns in a matrix */
  94.   void switch_rows(int row1, int row2);    /* swap two rows in a matrix */
  95.   void col_multiply_add(double multiplier, int reference_col, int changing_col);
  96.       /* changing_col += reference_col * multiplier */
  97.   void row_multiply_add(double multiplier, int reference_row, int changing_row);
  98.           /* changing_row += reference_row * multiplier */
  99.   double min();            /* find minimum element in the matrix */
  100.   double max();            /* find maximum element in the matrix */
  101.   double mean();        /* average all the elements of the matrix */
  102.   matrix & transpose();        /* transpose a square matrix */
  103.   double variance();        /* statistical variance of all elements */
  104.   matrix & scale();        /* Scale a matrix (used in L-U decomposition) */
  105.   void bitcopy(matrix& from, matrix& to); /* make an image of a matrix */
  106.   matrix & lu_decompose(matrix& indx, int& d );
  107.             /* Returns the L-U decomposition of a matrix */
  108.   void lu_back_subst(matrix& indx, matrix& b);
  109.             /* Uses lu-decomposition for matrix inverse */
  110.   void copy_column(matrix& m, int from_col, int to_col);
  111.   matrix & inverse();    /* Finally.  Invert a matrix. */
  112.   double determinant();
  113.   void error(char *s1) {  cerr << s1 << "\n";  exit(1); }
  114.   void error2(char *s1, char *s2 ) { cerr << s1 << " " << s2 << "\n"; exit(1); }
  115.   Cheb_vector & operator*(Cheb_vector & C);
  116. };
  117.  
  118.  
  119. /* 
  120. -*++ class RowVec: row vectors -- works with all sensible matrix operations
  121. ** 
  122. ** (*++ history: 
  123. **      8 Dec 87    Bruce Eckel    Creation date
  124. ** ++*)
  125. ** 
  126. ** (*++ detailed: 
  127. ** ++*)
  128. */
  129.  
  130. class RowVec: public matrix {
  131.   public:
  132.     RowVec() {;}
  133.     RowVec(int n) : (1,n,0) {;}
  134.     RowVec(int n, double v) : (1,n,v) {;} // initialize to v
  135.     double  & operator[](int x) { return (*(matrix *)this)[0][x]; }
  136.     int size() { return cols(); }
  137. };
  138.  
  139.  
  140.  
  141. /* 
  142. -*++ class ColVec: column vectors -- works with all sensible matrix operations
  143. ** 
  144. ** (*++ history: 
  145. **      8 Dec 87    Bruce Eckel    Creation date
  146. ** ++*)
  147. ** 
  148. ** (*++ detailed: 
  149. ** ++*)
  150. */
  151.  
  152. class ColVec: public matrix { 
  153.   public:
  154.     ColVec()  {;}
  155.     ColVec(int n) : (n,1,0) {;}
  156.     ColVec(int n, double v) : (n,1,v) {;} /* initialize to v */
  157.     ColVec(ColVec & x) : (x) {;}
  158.     double  & operator[](int x) { return (*(matrix *)this)[x][0]; }
  159.     double * doublvect();
  160.     int size() { return rows(); }
  161.     ColVec & operator=(ColVec & rval) 
  162.     { (matrix &)(*this) = (matrix &)rval; return (*this);}
  163. };
  164.  
  165.  
  166. /* 
  167. -*++ class indep_vec: A vector for holding the independent variable
  168. ** 
  169. ** (*++ history: 
  170. **     15 Dec 87    Bruce Eckel    Creation date
  171. ** ++*)
  172. ** 
  173. ** (*++ detailed: If the data is evenly spaced, it is stored as an offset and
  174. ** and increment.  If it is randomly spaced, it is stored directly.
  175. ** ++*)
  176. */
  177.  
  178. class indep_vec {
  179.      int ev_spc;        /* flag: independent variable is evenly spaced */
  180.      double strt_val;        
  181.      double delta;        
  182.      int Vsize;
  183.      ColVec *RandData;        /* if random, vector to store actual data */
  184.    public:
  185.      indep_vec(int sz, double v = 0) {        /* defaults to random */
  186.      RandData = new ColVec(sz,v);
  187.      Vsize = sz;
  188.      ev_spc = 0;
  189.      strt_val = delta = 0;
  190.      }
  191.      indep_vec(int sz, double start, double inc) {
  192.      ev_spc = 1;
  193.      Vsize = sz;
  194.      strt_val = start;
  195.      delta = inc;
  196.      }
  197.      ~indep_vec() { if(ev_spc) delete RandData; }
  198.      void IndepSet(int i, double val) { (*RandData)[i] = val;}
  199.      double Indep(int i) {return ev_spc? strt_val + i*delta : (*RandData)[i];}
  200.      int even_spaced() { return ev_spc; } /* info about this vector */
  201.      int random() { return !ev_spc; }
  202.      double offset() { return strt_val; }
  203.      double increment() { return delta; }
  204.      int size() { return Vsize; }
  205.  };
  206.  
  207. /* 
  208. -*++ class XY_vec: a column vector of associated X-Y data.
  209. ** 
  210. ** (*++ history: 
  211. **     10 Dec 87    Bruce Eckel    Creation date
  212. ** ++*)
  213. ** 
  214. ** (*++ detailed: Two associated vectors of data.  This is generally used in
  215. ** making x-y plots.  X is the independent variable.
  216. ** ++*)
  217. */
  218.  
  219.  class XY_vec {
  220.    public:
  221.      indep_vec X;
  222.      ColVec Y;
  223.      XY_vec(int n) : X(n,0), Y(n,0) {;}
  224.      XY_vec(int n, double v) : X(n,v), Y(n,v) {;} // initialize to v 
  225.      XY_vec(int n, double start, double inc) : X(n,start,inc), Y(n,0) {;}
  226.      int size() { return Y.rows(); } // length of XY vector
  227.      friend ostream& operator<<(ostream &s, XY_vec& v);
  228.      float * floatvect();
  229. //     XY_vec & operator+(XY_vec &);
  230. };
  231.