home *** CD-ROM | disk | FTP | other *** search
- #include <stream.h>
- #define TRUE 1
- #define FALSE 0
- #define EOL 10 /* end of line */
- #define TINY ((double)1e-20)
- #define NL "\n"
-
- /*
- -*++ class _Vector: thing to build matrices out of
- **
- ** (*++ history:
- ** 4 Dec 87 Bruce Eckel Creation date
- ** ++*)
- **
- ** (*++ detailed:
- ** ++*)
- */
-
- class _Vector {
- double *v_d;
- int v_sz;
- public:
- _Vector(int sz) { v_d = new double[v_sz = sz]; }
- ~_Vector() { delete v_d; }
- void error(int x) { cout << "vector index " << x << " out of bounds\n"; }
- double& operator[](int x) {
- if (x >= 0 && x < v_sz)
- return v_d[x];
- else { error(x); return v_d[x];}
- }
- };
-
- class Cheb_vector;
-
- /*
- -*++ class matrix: double matrix class. Rectangular, any size.
- **
- ** (*++ history:
- ** 4 Dec 87 Bruce Eckel Creation date. Was actually created
- ** long before this, but the GNUemacs documentation system was just made.
- ** ++*)
- **
- ** (*++ detailed: This class can be used by itself or can be derived from. It
- ** includes determinants, inversion, and several statistical functions as well
- ** as the usual matrix stuff. You can index into a matrix using the m[][]
- ** construct; this does bounds checking for you.
- ** ++*)
- */
-
- class matrix {
- struct matrep {
- _Vector** m_vec;
- int m_rows;
- int m_cols;
- int n; /* reference count */
- };
- matrep *p;
- public:
- matrix(int rows = 1, int columns = 1, double initval = 0);
- matrix(matrix& x);
- matrix(char * initfile); /* initialize from a "standard" matrix file */
- matrix(char * flag, int dimension); /* create an identity or random matrix */
- ~matrix();
- int rows() { return p->m_rows; }; /* max rows in matrix */
- int cols() { return p->m_cols; }; /* max cols in matrix */
- _Vector& operator[](int x) {
- if (x >= 0 && x < rows())
- return *(p->m_vec[x]);
- else {
- cerr << "row index out of bounds: " << x << "\n";
- return *(p->m_vec[x]);
- }
- }
- void data() { cout << "rows = " << rows() << " cols = " << cols()
- << " ref count = " << p->n << "\n"; };
- friend ostream& operator<<(ostream &s, matrix& m);
- friend istream& operator>>(istream &s, matrix& m);
- void write_standard(char * filename, char * msg = "");
- /* write matrix to `standard' file */
- void file_get(char * infile); /* read a matrix in from a file */
- void file_put(char * outfile); /* write a matrix out to a file */
- matrix & operator+(matrix&); /* add two matrices */
- matrix & operator+(double arg); /* add a scalar to each el */
- matrix & operator-(matrix& arg); /* subtract two matrices */
- matrix & operator-(); /* unary minus */
- matrix & operator*(double arg); /* multiply by a scalar */
- matrix & operator*(matrix& arg); /* matrix multiply */
- matrix & operator=(matrix&); /* matrix assignment */
- int operator==(matrix& arg); /* test for equivalence */
- int operator!=(matrix& arg); /* test for inequivalence */
- void set_row(int row, double value); /* set a row to a value */
- void set_column(int column, double value); /* set a column to a value */
- void switch_columns(int col1, int col2); /* swap two columns in a matrix */
- void switch_rows(int row1, int row2); /* swap two rows in a matrix */
- void col_multiply_add(double multiplier, int reference_col, int changing_col);
- /* changing_col += reference_col * multiplier */
- void row_multiply_add(double multiplier, int reference_row, int changing_row);
- /* changing_row += reference_row * multiplier */
- double min(); /* find minimum element in the matrix */
- double max(); /* find maximum element in the matrix */
- double mean(); /* average all the elements of the matrix */
- matrix & transpose(); /* transpose a square matrix */
- double variance(); /* statistical variance of all elements */
- matrix & scale(); /* Scale a matrix (used in L-U decomposition) */
- void bitcopy(matrix& from, matrix& to); /* make an image of a matrix */
- matrix & lu_decompose(matrix& indx, int& d );
- /* Returns the L-U decomposition of a matrix */
- void lu_back_subst(matrix& indx, matrix& b);
- /* Uses lu-decomposition for matrix inverse */
- void copy_column(matrix& m, int from_col, int to_col);
- matrix & inverse(); /* Finally. Invert a matrix. */
- double determinant();
- void error(char *s1) { cerr << s1 << "\n"; exit(1); }
- void error2(char *s1, char *s2 ) { cerr << s1 << " " << s2 << "\n"; exit(1); }
- Cheb_vector & operator*(Cheb_vector & C);
- };
-
-
- /*
- -*++ class RowVec: row vectors -- works with all sensible matrix operations
- **
- ** (*++ history:
- ** 8 Dec 87 Bruce Eckel Creation date
- ** ++*)
- **
- ** (*++ detailed:
- ** ++*)
- */
-
- class RowVec: public matrix {
- public:
- RowVec() {;}
- RowVec(int n) : (1,n,0) {;}
- RowVec(int n, double v) : (1,n,v) {;} // initialize to v
- double & operator[](int x) { return (*(matrix *)this)[0][x]; }
- int size() { return cols(); }
- };
-
-
-
- /*
- -*++ class ColVec: column vectors -- works with all sensible matrix operations
- **
- ** (*++ history:
- ** 8 Dec 87 Bruce Eckel Creation date
- ** ++*)
- **
- ** (*++ detailed:
- ** ++*)
- */
-
- class ColVec: public matrix {
- public:
- ColVec() {;}
- ColVec(int n) : (n,1,0) {;}
- ColVec(int n, double v) : (n,1,v) {;} /* initialize to v */
- ColVec(ColVec & x) : (x) {;}
- double & operator[](int x) { return (*(matrix *)this)[x][0]; }
- double * doublvect();
- int size() { return rows(); }
- ColVec & operator=(ColVec & rval)
- { (matrix &)(*this) = (matrix &)rval; return (*this);}
- };
-
-
- /*
- -*++ class indep_vec: A vector for holding the independent variable
- **
- ** (*++ history:
- ** 15 Dec 87 Bruce Eckel Creation date
- ** ++*)
- **
- ** (*++ detailed: If the data is evenly spaced, it is stored as an offset and
- ** and increment. If it is randomly spaced, it is stored directly.
- ** ++*)
- */
-
- class indep_vec {
- int ev_spc; /* flag: independent variable is evenly spaced */
- double strt_val;
- double delta;
- int Vsize;
- ColVec *RandData; /* if random, vector to store actual data */
- public:
- indep_vec(int sz, double v = 0) { /* defaults to random */
- RandData = new ColVec(sz,v);
- Vsize = sz;
- ev_spc = 0;
- strt_val = delta = 0;
- }
- indep_vec(int sz, double start, double inc) {
- ev_spc = 1;
- Vsize = sz;
- strt_val = start;
- delta = inc;
- }
- ~indep_vec() { if(ev_spc) delete RandData; }
- void IndepSet(int i, double val) { (*RandData)[i] = val;}
- double Indep(int i) {return ev_spc? strt_val + i*delta : (*RandData)[i];}
- int even_spaced() { return ev_spc; } /* info about this vector */
- int random() { return !ev_spc; }
- double offset() { return strt_val; }
- double increment() { return delta; }
- int size() { return Vsize; }
- };
-
- /*
- -*++ class XY_vec: a column vector of associated X-Y data.
- **
- ** (*++ history:
- ** 10 Dec 87 Bruce Eckel Creation date
- ** ++*)
- **
- ** (*++ detailed: Two associated vectors of data. This is generally used in
- ** making x-y plots. X is the independent variable.
- ** ++*)
- */
-
- class XY_vec {
- public:
- indep_vec X;
- ColVec Y;
- XY_vec(int n) : X(n,0), Y(n,0) {;}
- XY_vec(int n, double v) : X(n,v), Y(n,v) {;} // initialize to v
- XY_vec(int n, double start, double inc) : X(n,start,inc), Y(n,0) {;}
- int size() { return Y.rows(); } // length of XY vector
- friend ostream& operator<<(ostream &s, XY_vec& v);
- float * floatvect();
- // XY_vec & operator+(XY_vec &);
- };