Mercurial > hg > octave-nkf
changeset 2828:92826d6e8bd9
[project @ 1997-03-25 23:41:41 by jwe]
author | jwe |
---|---|
date | Tue, 25 Mar 1997 23:50:08 +0000 |
parents | 2387b8694c75 |
children | 6655bdca97fb |
files | liboctave/CDiagMatrix.h liboctave/CMatrix.cc liboctave/CMatrix.h liboctave/boolMatrix.cc liboctave/boolMatrix.h liboctave/chMatrix.cc liboctave/chMatrix.h liboctave/dMatrix.cc liboctave/dMatrix.h liboctave/idx-vector.cc liboctave/idx-vector.h liboctave/mx-base.h liboctave/mx-defs.h liboctave/mx-inlines.cc |
diffstat | 14 files changed, 375 insertions(+), 1667 deletions(-) [+] |
line wrap: on
line diff
--- a/liboctave/CDiagMatrix.h +++ b/liboctave/CDiagMatrix.h @@ -125,26 +125,12 @@ friend ComplexDiagMatrix operator * (const ComplexDiagMatrix& a, const ComplexDiagMatrix& b); - friend ComplexDiagMatrix operator + (const ComplexDiagMatrix& a, - const DiagMatrix& b); - friend ComplexDiagMatrix operator - (const ComplexDiagMatrix& a, - const DiagMatrix& b); friend ComplexDiagMatrix operator * (const ComplexDiagMatrix& a, const DiagMatrix& b); - friend ComplexDiagMatrix operator + (const DiagMatrix& a, - const ComplexDiagMatrix& b); - friend ComplexDiagMatrix operator - (const DiagMatrix& a, - const ComplexDiagMatrix& b); friend ComplexDiagMatrix operator * (const DiagMatrix& a, const ComplexDiagMatrix& b); - friend ComplexDiagMatrix product (const ComplexDiagMatrix& a, - const DiagMatrix& b); - - friend ComplexDiagMatrix product (const DiagMatrix& a, - const ComplexDiagMatrix& b); - // other operations ComplexColumnVector diag (void) const;
--- a/liboctave/CMatrix.cc +++ b/liboctave/CMatrix.cc @@ -38,6 +38,7 @@ #include <sys/types.h> #endif +#include "CMatrix.h" #include "CmplxAEPBAL.h" #include "CmplxDET.h" #include "CmplxSCHUR.h" @@ -48,6 +49,8 @@ #include "lo-mappers.h" #include "lo-utils.h" #include "mx-base.h" +#include "mx-cm-dm.h" +#include "mx-cm-s.h" #include "mx-inlines.cc" #include "oct-cmplx.h" @@ -159,6 +162,13 @@ // XXX FIXME XXX -- could we use a templated mixed-type copy function // here? +ComplexMatrix::ComplexMatrix (const boolMatrix& a) +{ + for (int i = 0; i < a.cols (); i++) + for (int j = 0; j < a.rows (); j++) + elem (i, j) = a.elem (i, j); +} + ComplexMatrix::ComplexMatrix (const charMatrix& a) { for (int i = 0; i < a.cols (); i++) @@ -1712,94 +1722,6 @@ return retval; } -// diagonal matrix by scalar -> matrix operations - -ComplexMatrix -operator + (const DiagMatrix& a, const Complex& s) -{ - ComplexMatrix tmp (a.rows (), a.cols (), s); - return a + tmp; -} - -ComplexMatrix -operator - (const DiagMatrix& a, const Complex& s) -{ - ComplexMatrix tmp (a.rows (), a.cols (), -s); - return a + tmp; -} - -ComplexMatrix -operator + (const ComplexDiagMatrix& a, double s) -{ - ComplexMatrix tmp (a.rows (), a.cols (), s); - return a + tmp; -} - -ComplexMatrix -operator - (const ComplexDiagMatrix& a, double s) -{ - ComplexMatrix tmp (a.rows (), a.cols (), -s); - return a + tmp; -} - -ComplexMatrix -operator + (const ComplexDiagMatrix& a, const Complex& s) -{ - ComplexMatrix tmp (a.rows (), a.cols (), s); - return a + tmp; -} - -ComplexMatrix -operator - (const ComplexDiagMatrix& a, const Complex& s) -{ - ComplexMatrix tmp (a.rows (), a.cols (), -s); - return a + tmp; -} - -// scalar by diagonal matrix -> matrix operations - -ComplexMatrix -operator + (const Complex& s, const DiagMatrix& a) -{ - ComplexMatrix tmp (a.rows (), a.cols (), s); - return tmp + a; -} - -ComplexMatrix -operator - (const Complex& s, const DiagMatrix& a) -{ - ComplexMatrix tmp (a.rows (), a.cols (), s); - return tmp - a; -} - -ComplexMatrix -operator + (double s, const ComplexDiagMatrix& a) -{ - ComplexMatrix tmp (a.rows (), a.cols (), s); - return tmp + a; -} - -ComplexMatrix -operator - (double s, const ComplexDiagMatrix& a) -{ - ComplexMatrix tmp (a.rows (), a.cols (), s); - return tmp - a; -} - -ComplexMatrix -operator + (const Complex& s, const ComplexDiagMatrix& a) -{ - ComplexMatrix tmp (a.rows (), a.cols (), s); - return tmp + a; -} - -ComplexMatrix -operator - (const Complex& s, const ComplexDiagMatrix& a) -{ - ComplexMatrix tmp (a.rows (), a.cols (), s); - return tmp - a; -} - // matrix by diagonal matrix -> matrix operations ComplexMatrix& @@ -1886,411 +1808,6 @@ return *this; } -ComplexMatrix -operator + (const Matrix& m, const ComplexDiagMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator +", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (m); - for (int i = 0; i < a.length (); i++) - result.elem (i, i) += a.elem (i, i); - - return result; -} - -ComplexMatrix -operator - (const Matrix& m, const ComplexDiagMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator -", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (m); - for (int i = 0; i < a.length (); i++) - result.elem (i, i) -= a.elem (i, i); - - return result; -} - -ComplexMatrix -operator * (const Matrix& m, const ComplexDiagMatrix& a) -{ - ComplexMatrix retval; - - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nc != a_nr) - gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); - else - { - if (nr == 0 || nc == 0 || a_nc == 0) - retval.resize (nr, a_nc, 0.0); - else - { - retval.resize (nr, a_nc); - Complex *c = retval.fortran_vec (); - - Complex *ctmp = 0; - - for (int j = 0; j < a.length (); j++) - { - int idx = j * nr; - ctmp = c + idx; - if (a.elem (j, j) == 1.0) - { - for (int i = 0; i < nr; i++) - ctmp[i] = m.elem (i, j); - } - else if (a.elem (j, j) == 0.0) - { - for (int i = 0; i < nr; i++) - ctmp[i] = 0.0; - } - else - { - for (int i = 0; i < nr; i++) - ctmp[i] = a.elem (j, j) * m.elem (i, j); - } - } - - if (a_nr < a_nc) - { - for (int i = nr * nc; i < nr * a_nc; i++) - ctmp[i] = 0.0; - } - } - } - - return retval; -} - -// diagonal matrix by matrix -> matrix operations - -ComplexMatrix -operator + (const DiagMatrix& m, const ComplexMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator +", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (a); - for (int i = 0; i < m.length (); i++) - result.elem (i, i) += m.elem (i, i); - - return result; -} - -ComplexMatrix -operator - (const DiagMatrix& m, const ComplexMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator -", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (-a); - for (int i = 0; i < m.length (); i++) - result.elem (i, i) += m.elem (i, i); - - return result; -} - -ComplexMatrix -operator * (const DiagMatrix& m, const ComplexMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nc != a_nr) - { - gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0 || a_nc == 0) - return ComplexMatrix (nr, nc, 0.0); - - ComplexMatrix c (nr, a_nc); - - for (int i = 0; i < m.length (); i++) - { - if (m.elem (i, i) == 1.0) - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = a.elem (i, j); - } - else if (m.elem (i, i) == 0.0) - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = 0.0; - } - else - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = m.elem (i, i) * a.elem (i, j); - } - } - - if (nr > nc) - { - for (int j = 0; j < a_nc; j++) - for (int i = a_nr; i < nr; i++) - c.elem (i, j) = 0.0; - } - - return c; -} - -ComplexMatrix -operator + (const ComplexDiagMatrix& m, const Matrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator +", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (a); - for (int i = 0; i < m.length (); i++) - result.elem (i, i) += m.elem (i, i); - - return result; -} - -ComplexMatrix -operator - (const ComplexDiagMatrix& m, const Matrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator -", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (-a); - for (int i = 0; i < m.length (); i++) - result.elem (i, i) += m.elem (i, i); - - return result; -} - -ComplexMatrix -operator * (const ComplexDiagMatrix& m, const Matrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nc != a_nr) - { - gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0 || a_nc == 0) - return ComplexMatrix (nr, a_nc, 0.0); - - ComplexMatrix c (nr, a_nc); - - for (int i = 0; i < m.length (); i++) - { - if (m.elem (i, i) == 1.0) - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = a.elem (i, j); - } - else if (m.elem (i, i) == 0.0) - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = 0.0; - } - else - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = m.elem (i, i) * a.elem (i, j); - } - } - - if (nr > nc) - { - for (int j = 0; j < a_nc; j++) - for (int i = a_nr; i < nr; i++) - c.elem (i, j) = 0.0; - } - - return c; -} - -ComplexMatrix -operator + (const ComplexDiagMatrix& m, const ComplexMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator +", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (a); - for (int i = 0; i < m.length (); i++) - result.elem (i, i) += m.elem (i, i); - - return result; -} - -ComplexMatrix -operator - (const ComplexDiagMatrix& m, const ComplexMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator -", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (-a); - for (int i = 0; i < m.length (); i++) - result.elem (i, i) += m.elem (i, i); - - return result; -} - -ComplexMatrix -operator * (const ComplexDiagMatrix& m, const ComplexMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nc != a_nr) - { - gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0 || a_nc == 0) - return ComplexMatrix (nr, a_nc, 0.0); - - ComplexMatrix c (nr, a_nc); - - for (int i = 0; i < m.length (); i++) - { - if (m.elem (i, i) == 1.0) - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = a.elem (i, j); - } - else if (m.elem (i, i) == 0.0) - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = 0.0; - } - else - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = m.elem (i, i) * a.elem (i, j); - } - } - - if (nr > nc) - { - for (int j = 0; j < a_nc; j++) - for (int i = a_nr; i < nr; i++) - c.elem (i, j) = 0.0; - } - - return c; -} - // matrix by matrix -> matrix operations ComplexMatrix& @@ -2397,553 +1914,6 @@ return Matrix (not (data (), length ()), rows (), cols ()); } -// matrix by scalar -> matrix operations - -ComplexMatrix -operator + (const Matrix& a, const Complex& s) -{ - return ComplexMatrix (add (a.data (), a.length (), s), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator - (const Matrix& a, const Complex& s) -{ - return ComplexMatrix (subtract (a.data (), a.length (), s), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator * (const Matrix& a, const Complex& s) -{ - return ComplexMatrix (multiply (a.data (), a.length (), s), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator / (const Matrix& a, const Complex& s) -{ - return ComplexMatrix (divide (a.data (), a.length (), s), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator + (const ComplexMatrix& a, double s) -{ - return ComplexMatrix (add (a.data (), a.length (), s), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator - (const ComplexMatrix& a, double s) -{ - return ComplexMatrix (subtract (a.data (), a.length (), s), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator * (const ComplexMatrix& a, double s) -{ - return ComplexMatrix (multiply (a.data (), a.length (), s), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator / (const ComplexMatrix& a, double s) -{ - return ComplexMatrix (divide (a.data (), a.length (), s), - a.rows (), a.cols ()); -} - -// scalar by matrix -> matrix operations - -ComplexMatrix -operator + (double s, const ComplexMatrix& a) -{ - return ComplexMatrix (add (a.data (), a.length (), s), a.rows (), - a.cols ()); -} - -ComplexMatrix -operator - (double s, const ComplexMatrix& a) -{ - return ComplexMatrix (subtract (s, a.data (), a.length ()), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator * (double s, const ComplexMatrix& a) -{ - return ComplexMatrix (multiply (a.data (), a.length (), s), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator / (double s, const ComplexMatrix& a) -{ - return ComplexMatrix (divide (s, a.data (), a.length ()), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator + (const Complex& s, const Matrix& a) -{ - return ComplexMatrix (add (s, a.data (), a.length ()), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator - (const Complex& s, const Matrix& a) -{ - return ComplexMatrix (subtract (s, a.data (), a.length ()), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator * (const Complex& s, const Matrix& a) -{ - return ComplexMatrix (multiply (a.data (), a.length (), s), - a.rows (), a.cols ()); -} - -ComplexMatrix -operator / (const Complex& s, const Matrix& a) -{ - return ComplexMatrix (divide (s, a.data (), a.length ()), - a.rows (), a.cols ()); -} - -// matrix by diagonal matrix -> matrix operations - -ComplexMatrix -operator + (const ComplexMatrix& m, const DiagMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator +", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (m); - for (int i = 0; i < a.length (); i++) - result.elem (i, i) += a.elem (i, i); - - return result; -} - -ComplexMatrix -operator - (const ComplexMatrix& m, const DiagMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator -", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (m); - for (int i = 0; i < a.length (); i++) - result.elem (i, i) -= a.elem (i, i); - - return result; -} - -ComplexMatrix -operator * (const ComplexMatrix& m, const DiagMatrix& a) -{ - ComplexMatrix retval; - - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nc != a_nr) - gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); - else - { - if (nr == 0 || nc == 0 || a_nc == 0) - retval.resize (nr, nc, 0.0); - else - { - retval.resize (nr, a_nc); - Complex *c = retval.fortran_vec (); - Complex *ctmp = 0; - - for (int j = 0; j < a.length (); j++) - { - int idx = j * nr; - ctmp = c + idx; - if (a.elem (j, j) == 1.0) - { - for (int i = 0; i < nr; i++) - ctmp[i] = m.elem (i, j); - } - else if (a.elem (j, j) == 0.0) - { - for (int i = 0; i < nr; i++) - ctmp[i] = 0.0; - } - else - { - for (int i = 0; i < nr; i++) - ctmp[i] = a.elem (j, j) * m.elem (i, j); - } - } - - if (a.rows () < a_nc) - { - for (int i = nr * nc; i < nr * a_nc; i++) - ctmp[i] = 0.0; - } - } - } - - return retval; -} - -ComplexMatrix -operator + (const ComplexMatrix& m, const ComplexDiagMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator +", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (m); - for (int i = 0; i < a.length (); i++) - result.elem (i, i) += a.elem (i, i); - - return result; -} - -ComplexMatrix -operator - (const ComplexMatrix& m, const ComplexDiagMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator -", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - ComplexMatrix result (m); - for (int i = 0; i < a.length (); i++) - result.elem (i, i) -= a.elem (i, i); - - return result; -} - -ComplexMatrix -operator * (const ComplexMatrix& m, const ComplexDiagMatrix& a) -{ - ComplexMatrix retval; - - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nc != a_nr) - gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); - else - { - if (nr == 0 || nc == 0 || a_nc == 0) - retval.resize (nr, nc, 0.0); - else - { - retval.resize (nr, nc); - Complex *c = retval.fortran_vec (); - Complex *ctmp = 0; - - for (int j = 0; j < a.length (); j++) - { - int idx = j * nr; - ctmp = c + idx; - if (a.elem (j, j) == 1.0) - { - for (int i = 0; i < nr; i++) - ctmp[i] = m.elem (i, j); - } - else if (a.elem (j, j) == 0.0) - { - for (int i = 0; i < nr; i++) - ctmp[i] = 0.0; - } - else - { - for (int i = 0; i < nr; i++) - ctmp[i] = a.elem (j, j) * m.elem (i, j); - } - } - - if (a.rows () < a_nc) - { - for (int i = nr * nc; i < nr * a_nc; i++) - ctmp[i] = 0.0; - } - } - } - - return retval; -} - -// matrix by matrix -> matrix operations - -ComplexMatrix -operator + (const ComplexMatrix& m, const Matrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator +", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - return ComplexMatrix (add (m.data (), a.data (), m.length ()), nr, nc); -} - -ComplexMatrix -operator - (const ComplexMatrix& m, const Matrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator -", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - return ComplexMatrix (subtract (m.data (), a.data (), m.length ()), nr, nc); -} - -ComplexMatrix -operator + (const Matrix& m, const ComplexMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator +", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - return ComplexMatrix (add (m.data (), a.data (), m.length ()), nr, nc); -} - -ComplexMatrix -operator - (const Matrix& m, const ComplexMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator -", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - return ComplexMatrix (subtract (m.data (), a.data (), m.length ()), nr, nc); -} - -ComplexMatrix -operator * (const ComplexMatrix& m, const Matrix& a) -{ - ComplexMatrix tmp (a); - return m * tmp; -} - -ComplexMatrix -operator * (const Matrix& m, const ComplexMatrix& a) -{ - ComplexMatrix tmp (m); - return tmp * a; -} - -ComplexMatrix -operator * (const ComplexMatrix& m, const ComplexMatrix& a) -{ - ComplexMatrix retval; - - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nc != a_nr) - gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); - else - { - if (nr == 0 || nc == 0 || a_nc == 0) - retval.resize (nr, nc, 0.0); - else - { - int ld = nr; - int lda = a.rows (); - - retval.resize (nr, a_nc); - Complex *c = retval.fortran_vec (); - - F77_XFCN (zgemm, ZGEMM, ("N", "N", nr, a_nc, nc, 1.0, - m.data (), ld, a.data (), lda, 0.0, - c, nr, 1L, 1L)); - - if (f77_exception_encountered) - (*current_liboctave_error_handler) - ("unrecoverable error in zgemm"); - } - } - - return retval; -} - -ComplexMatrix -product (const ComplexMatrix& m, const Matrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("product", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - return ComplexMatrix (multiply (m.data (), a.data (), m.length ()), nr, nc); -} - -ComplexMatrix -quotient (const ComplexMatrix& m, const Matrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("quotient", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - return ComplexMatrix (divide (m.data (), a.data (), m.length ()), nr, nc); -} - -ComplexMatrix -product (const Matrix& m, const ComplexMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("product", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - return ComplexMatrix (multiply (m.data (), a.data (), m.length ()), nr, nc); -} - -ComplexMatrix -quotient (const Matrix& m, const ComplexMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("quotient", nr, nc, a_nr, a_nc); - return ComplexMatrix (); - } - - if (nr == 0 || nc == 0) - return ComplexMatrix (nr, nc); - - return ComplexMatrix (divide (m.data (), a.data (), m.length ()), nr, nc); -} - // other operations ComplexMatrix @@ -3856,6 +2826,58 @@ return retval; } +ComplexMatrix +operator * (const ComplexMatrix& m, const Matrix& a) +{ + ComplexMatrix tmp (a); + return m * tmp; +} + +ComplexMatrix +operator * (const Matrix& m, const ComplexMatrix& a) +{ + ComplexMatrix tmp (m); + return tmp * a; +} + +ComplexMatrix +operator * (const ComplexMatrix& m, const ComplexMatrix& a) +{ + ComplexMatrix retval; + + int nr = m.rows (); + int nc = m.cols (); + + int a_nr = a.rows (); + int a_nc = a.cols (); + + if (nc != a_nr) + gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); + else + { + if (nr == 0 || nc == 0 || a_nc == 0) + retval.resize (nr, nc, 0.0); + else + { + int ld = nr; + int lda = a.rows (); + + retval.resize (nr, a_nc); + Complex *c = retval.fortran_vec (); + + F77_XFCN (zgemm, ZGEMM, ("N", "N", nr, a_nc, nc, 1.0, + m.data (), ld, a.data (), lda, 0.0, + c, nr, 1L, 1L)); + + if (f77_exception_encountered) + (*current_liboctave_error_handler) + ("unrecoverable error in zgemm"); + } + } + + return retval; +} + /* ;;; Local Variables: *** ;;; mode: C++ ***
--- a/liboctave/CMatrix.h +++ b/liboctave/CMatrix.h @@ -61,6 +61,7 @@ ComplexMatrix (const ComplexColumnVector& cv); ComplexMatrix (const ComplexDiagMatrix& a); + ComplexMatrix (const boolMatrix& a); ComplexMatrix (const charMatrix& a); ComplexMatrix& operator = (const ComplexMatrix& a) @@ -180,32 +181,6 @@ friend ComplexMatrix operator * (const ComplexColumnVector& a, const ComplexRowVector& b); - // diagonal matrix by scalar -> matrix operations - - friend ComplexMatrix operator + (const DiagMatrix& a, const Complex& s); - friend ComplexMatrix operator - (const DiagMatrix& a, const Complex& s); - - friend ComplexMatrix operator + (const ComplexDiagMatrix& a, double s); - friend ComplexMatrix operator - (const ComplexDiagMatrix& a, double s); - - friend ComplexMatrix operator + (const ComplexDiagMatrix& a, - const Complex& s); - friend ComplexMatrix operator - (const ComplexDiagMatrix& a, - const Complex& s); - - // scalar by diagonal matrix -> matrix operations - - friend ComplexMatrix operator + (const Complex& s, const DiagMatrix& a); - friend ComplexMatrix operator - (const Complex& s, const DiagMatrix& a); - - friend ComplexMatrix operator + (double s, const ComplexDiagMatrix& a); - friend ComplexMatrix operator - (double s, const ComplexDiagMatrix& a); - - friend ComplexMatrix operator + (const Complex& s, - const ComplexDiagMatrix& a); - friend ComplexMatrix operator - (const Complex& s, - const ComplexDiagMatrix& a); - // matrix by diagonal matrix -> matrix operations ComplexMatrix& operator += (const DiagMatrix& a); @@ -214,36 +189,6 @@ ComplexMatrix& operator += (const ComplexDiagMatrix& a); ComplexMatrix& operator -= (const ComplexDiagMatrix& a); - friend ComplexMatrix operator + (const Matrix& a, - const ComplexDiagMatrix& b); - friend ComplexMatrix operator - (const Matrix& a, - const ComplexDiagMatrix& b); - friend ComplexMatrix operator * (const Matrix& a, - const ComplexDiagMatrix& b); - - // diagonal matrix by matrix -> matrix operations - - friend ComplexMatrix operator + (const DiagMatrix& a, - const ComplexMatrix& b); - friend ComplexMatrix operator - (const DiagMatrix& a, - const ComplexMatrix& b); - friend ComplexMatrix operator * (const DiagMatrix& a, - const ComplexMatrix& b); - - friend ComplexMatrix operator + (const ComplexDiagMatrix& a, - const Matrix& b); - friend ComplexMatrix operator - (const ComplexDiagMatrix& a, - const Matrix& b); - friend ComplexMatrix operator * (const ComplexDiagMatrix& a, - const Matrix& b); - - friend ComplexMatrix operator + (const ComplexDiagMatrix& a, - const ComplexMatrix& b); - friend ComplexMatrix operator - (const ComplexDiagMatrix& a, - const ComplexMatrix& b); - friend ComplexMatrix operator * (const ComplexDiagMatrix& a, - const ComplexMatrix& b); - // matrix by matrix -> matrix operations ComplexMatrix& operator += (const Matrix& a); @@ -256,67 +201,6 @@ Matrix operator ! (void) const; - // matrix by scalar -> matrix operations - - friend ComplexMatrix operator + (const Matrix& a, const Complex& s); - friend ComplexMatrix operator - (const Matrix& a, const Complex& s); - friend ComplexMatrix operator * (const Matrix& a, const Complex& s); - friend ComplexMatrix operator / (const Matrix& a, const Complex& s); - - friend ComplexMatrix operator + (const ComplexMatrix& a, double s); - friend ComplexMatrix operator - (const ComplexMatrix& a, double s); - friend ComplexMatrix operator * (const ComplexMatrix& a, double s); - friend ComplexMatrix operator / (const ComplexMatrix& a, double s); - - // scalar by matrix -> matrix operations - - friend ComplexMatrix operator + (double s, const ComplexMatrix& a); - friend ComplexMatrix operator - (double s, const ComplexMatrix& a); - friend ComplexMatrix operator * (double s, const ComplexMatrix& a); - friend ComplexMatrix operator / (double s, const ComplexMatrix& a); - - friend ComplexMatrix operator + (const Complex& s, const Matrix& a); - friend ComplexMatrix operator - (const Complex& s, const Matrix& a); - friend ComplexMatrix operator * (const Complex& s, const Matrix& a); - friend ComplexMatrix operator / (const Complex& s, const Matrix& a); - - // matrix by diagonal matrix -> matrix operations - - friend ComplexMatrix operator + (const ComplexMatrix& a, - const DiagMatrix& b); - friend ComplexMatrix operator - (const ComplexMatrix& a, - const DiagMatrix& b); - friend ComplexMatrix operator * (const ComplexMatrix& a, - const DiagMatrix& b); - - friend ComplexMatrix operator + (const ComplexMatrix& a, - const ComplexDiagMatrix& b); - friend ComplexMatrix operator - (const ComplexMatrix& a, - const ComplexDiagMatrix& b); - friend ComplexMatrix operator * (const ComplexMatrix& a, - const ComplexDiagMatrix& b); - - // matrix by matrix -> matrix operations - - friend ComplexMatrix operator + (const ComplexMatrix& a, const Matrix& b); - friend ComplexMatrix operator - (const ComplexMatrix& a, const Matrix& b); - - friend ComplexMatrix operator + (const Matrix& a, const ComplexMatrix& b); - friend ComplexMatrix operator - (const Matrix& a, const ComplexMatrix& b); - - friend ComplexMatrix operator * (const ComplexMatrix& a, const Matrix& b); - - friend ComplexMatrix operator * (const Matrix& a, const ComplexMatrix& b); - - friend ComplexMatrix operator * (const ComplexMatrix& a, - const ComplexMatrix& b); - - friend ComplexMatrix product (const ComplexMatrix& a, const Matrix& b); - friend ComplexMatrix quotient (const ComplexMatrix& a, const Matrix& b); - - friend ComplexMatrix product (const Matrix& a, const ComplexMatrix& b); - friend ComplexMatrix quotient (const Matrix& a, const ComplexMatrix& b); - // other operations ComplexMatrix map (c_c_Mapper f) const; @@ -371,6 +255,10 @@ ComplexMatrix Sylvester (const ComplexMatrix&, const ComplexMatrix&, const ComplexMatrix&); +extern ComplexMatrix operator * (const Matrix&, const ComplexMatrix&); +extern ComplexMatrix operator * (const ComplexMatrix&, const Matrix&); +extern ComplexMatrix operator * (const ComplexMatrix&, const ComplexMatrix&); + #endif /*
new file mode 100644 --- /dev/null +++ b/liboctave/boolMatrix.cc @@ -0,0 +1,82 @@ +// Matrix manipulations. +/* + +Copyright (C) 1996 John W. Eaton + +This file is part of Octave. + +Octave is free software; you can redistribute it and/or modify it +under the terms of the GNU General Public License as published by the +Free Software Foundation; either version 2, or (at your option) any +later version. + +Octave is distributed in the hope that it will be useful, but WITHOUT +ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +for more details. + +You should have received a copy of the GNU General Public License +along with Octave; see the file COPYING. If not, write to the Free +Software Foundation, 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. + +*/ + +#if defined (__GNUG__) +#pragma implementation +#endif + +#ifdef HAVE_CONFIG_H +#include <config.h> +#endif + +#include <iostream.h> + +#include "lo-error.h" +#include "str-vec.h" +#include "mx-base.h" +#include "mx-inlines.cc" + +// boolMatrix class. + +bool +boolMatrix::operator == (const boolMatrix& a) const +{ + if (rows () != a.rows () || cols () != a.cols ()) + return 0; + + return equal (data (), a.data (), length ()); +} + +bool +boolMatrix::operator != (const boolMatrix& a) const +{ + return !(*this == a); +} + +boolMatrix& +boolMatrix::insert (const boolMatrix& a, int r, int c) +{ + Array2<bool>::insert (a, r, c); + return *this; +} + +boolMatrix +boolMatrix::transpose (void) const +{ + int nr = rows (); + int nc = cols (); + boolMatrix result (nc, nr); + if (length () > 0) + { + for (int j = 0; j < nc; j++) + for (int i = 0; i < nr; i++) + result.elem (j, i) = elem (i, j); + } + return result; +} + +/* +;;; Local Variables: *** +;;; mode: C++ *** +;;; End: *** +*/
new file mode 100644 --- /dev/null +++ b/liboctave/boolMatrix.h @@ -0,0 +1,78 @@ +/* + +Copyright (C) 1996 John W. Eaton + +This file is part of Octave. + +Octave is free software; you can redistribute it and/or modify it +under the terms of the GNU General Public License as published by the +Free Software Foundation; either version 2, or (at your option) any +later version. + +Octave is distributed in the hope that it will be useful, but WITHOUT +ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +for more details. + +You should have received a copy of the GNU General Public License +along with Octave; see the file COPYING. If not, write to the Free +Software Foundation, 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. + +*/ + +#if !defined (octave_boolMatrix_int_h) +#define octave_boolMatrix_int_h 1 + +#if defined (__GNUG__) +#pragma interface +#endif + +#include "MArray2.h" + +#include "mx-defs.h" + +class +boolMatrix : public MArray2<bool> +{ +public: + + boolMatrix (void) : MArray2<bool> () { } + boolMatrix (int r, int c) : MArray2<bool> (r, c) { } + boolMatrix (int r, int c, bool val) : MArray2<bool> (r, c, val) { } + boolMatrix (const MArray2<bool>& a) : MArray2<bool> (a) { } + boolMatrix (const boolMatrix& a) : MArray2<bool> (a) { } + + boolMatrix& operator = (const boolMatrix& a) + { + MArray2<bool>::operator = (a); + return *this; + } + + bool operator == (const boolMatrix& a) const; + bool operator != (const boolMatrix& a) const; + + // destructive insert/delete/reorder operations + + boolMatrix& insert (const boolMatrix& a, int r, int c); + + boolMatrix transpose (void) const; + +#if 0 + // i/o + + friend ostream& operator << (ostream& os, const Matrix& a); + friend istream& operator >> (istream& is, Matrix& a); +#endif + +private: + + boolMatrix (bool *b, int r, int c) : MArray2<bool> (b, r, c) { } +}; + +#endif + +/* +;;; Local Variables: *** +;;; mode: C++ *** +;;; End: *** +*/
--- a/liboctave/chMatrix.cc +++ b/liboctave/chMatrix.cc @@ -29,15 +29,10 @@ #include <config.h> #endif -#include <cstdio> -#include <cstring> - #include <string> #include <iostream.h> -// #include <sys/types.h> // XXX FIXME XXX - #include "lo-error.h" #include "str-vec.h" #include "mx-base.h"
--- a/liboctave/chMatrix.h +++ b/liboctave/chMatrix.h @@ -27,9 +27,6 @@ #pragma interface #endif -// For FILE... -#include <cstdio> - #include <string> #include "MArray2.h" @@ -69,172 +66,13 @@ string row_as_string (int r, bool strip_trailing_whitespace = false) const; -#if 0 - Matrix& insert (const RowVector& a, int r, int c); - Matrix& insert (const ColumnVector& a, int r, int c); - Matrix& insert (const DiagMatrix& a, int r, int c); - - Matrix& fill (char val); - Matrix& fill (char val, int r1, int c1, int r2, int c2); - - Matrix append (const Matrix& a) const; - Matrix append (const RowVector& a) const; - Matrix append (const ColumnVector& a) const; - Matrix append (const DiagMatrix& a) const; - - Matrix stack (const Matrix& a) const; - Matrix stack (const RowVector& a) const; - Matrix stack (const ColumnVector& a) const; - Matrix stack (const DiagMatrix& a) const; -#endif - charMatrix transpose (void) const; #if 0 - friend Matrix real (const ComplexMatrix& a); - friend Matrix imag (const ComplexMatrix& a); - - // resize is the destructive equivalent for this one - - Matrix extract (int r1, int c1, int r2, int c2) const; - - // extract row or column i. - - RowVector row (int i) const; - RowVector row (char *s) const; - - ColumnVector column (int i) const; - ColumnVector column (char *s) const; - - Matrix inverse (void) const; - Matrix inverse (int& info) const; - Matrix inverse (int& info, double& rcond) const; - - Matrix pseudo_inverse (double tol = 0.0); - - ComplexMatrix fourier (void) const; - ComplexMatrix ifourier (void) const; - - ComplexMatrix fourier2d (void) const; - ComplexMatrix ifourier2d (void) const; - - DET determinant (void) const; - DET determinant (int& info) const; - DET determinant (int& info, double& rcond) const; - - Matrix solve (const Matrix& b) const; - Matrix solve (const Matrix& b, int& info) const; - Matrix solve (const Matrix& b, int& info, double& rcond) const; - - ComplexMatrix solve (const ComplexMatrix& b) const; - ComplexMatrix solve (const ComplexMatrix& b, int& info) const; - ComplexMatrix solve (const ComplexMatrix& b, int& info, double& rcond) const; - - ColumnVector solve (const ColumnVector& b) const; - ColumnVector solve (const ColumnVector& b, int& info) const; - ColumnVector solve (const ColumnVector& b, int& info, double& rcond) const; - - ComplexColumnVector solve (const ComplexColumnVector& b) const; - ComplexColumnVector solve (const ComplexColumnVector& b, int& info) const; - ComplexColumnVector solve (const ComplexColumnVector& b, int& info, - double& rcond) const; - - Matrix lssolve (const Matrix& b) const; - Matrix lssolve (const Matrix& b, int& info) const; - Matrix lssolve (const Matrix& b, int& info, int& rank) const; - - ComplexMatrix lssolve (const ComplexMatrix& b) const; - ComplexMatrix lssolve (const ComplexMatrix& b, int& info) const; - ComplexMatrix lssolve (const ComplexMatrix& b, int& info, - int& rank) const; - - ColumnVector lssolve (const ColumnVector& b) const; - ColumnVector lssolve (const ColumnVector& b, int& info) const; - ColumnVector lssolve (const ColumnVector& b, int& info, int& rank) const; - - ComplexColumnVector lssolve (const ComplexColumnVector& b) const; - ComplexColumnVector lssolve (const ComplexColumnVector& b, int& info) const; - ComplexColumnVector lssolve (const ComplexColumnVector& b, int& info, - int& rank) const; - - Matrix& operator += (const Matrix& a); - Matrix& operator -= (const Matrix& a); - - Matrix& operator += (const DiagMatrix& a); - Matrix& operator -= (const DiagMatrix& a); - - // unary operations - - Matrix operator ! (void) const; - - // column vector by row vector -> matrix operations - - friend Matrix operator * (const ColumnVector& a, const RowVector& a); - - // diagonal matrix by scalar -> matrix operations - - friend Matrix operator + (const DiagMatrix& a, double s); - friend Matrix operator - (const DiagMatrix& a, double s); - - // scalar by diagonal matrix -> matrix operations - - friend Matrix operator + (double s, const DiagMatrix& a); - friend Matrix operator - (double s, const DiagMatrix& a); - - // matrix by diagonal matrix -> matrix operations - - friend Matrix operator + (const Matrix& a, const DiagMatrix& b); - friend Matrix operator - (const Matrix& a, const DiagMatrix& b); - friend Matrix operator * (const Matrix& a, const DiagMatrix& b); - - // diagonal matrix by matrix -> matrix operations - - friend Matrix operator + (const DiagMatrix& a, const Matrix& b); - friend Matrix operator - (const DiagMatrix& a, const Matrix& b); - friend Matrix operator * (const DiagMatrix& a, const Matrix& b); - - // matrix by matrix -> matrix operations - - friend Matrix operator * (const Matrix& a, const Matrix& b); - - // other operations - - friend Matrix map (d_d_Mapper f, const Matrix& a); - friend Matrix map (d_c_Mapper f, const ComplexMatrix& a); - - void map (d_d_Mapper f); - - Matrix all (void) const; - Matrix any (void) const; - - Matrix cumprod (void) const; - Matrix cumsum (void) const; - Matrix prod (void) const; - Matrix sum (void) const; - Matrix sumsq (void) const; - - ColumnVector diag (void) const; - ColumnVector diag (int k) const; - - ColumnVector row_min (void) const; - ColumnVector row_min_loc (void) const; - - ColumnVector row_max (void) const; - ColumnVector row_max_loc (void) const; - - RowVector column_min (void) const; - RowVector column_min_loc (void) const; - - RowVector column_max (void) const; - RowVector column_max_loc (void) const; - // i/o friend ostream& operator << (ostream& os, const Matrix& a); friend istream& operator >> (istream& is, Matrix& a); - - int read (FILE *fptr, const char *type); - int write (FILE *fptr, const char *type); #endif private:
--- a/liboctave/dMatrix.cc +++ b/liboctave/dMatrix.cc @@ -34,6 +34,7 @@ #include <iostream.h> #include "byte-swap.h" +#include "dMatrix.h" #include "dbleAEPBAL.h" #include "dbleDET.h" #include "dbleSCHUR.h" @@ -44,6 +45,7 @@ #include "lo-mappers.h" #include "lo-utils.h" #include "mx-base.h" +#include "mx-m-dm.h" #include "mx-inlines.cc" #include "oct-cmplx.h" @@ -135,6 +137,14 @@ // XXX FIXME XXX -- could we use a templated mixed-type copy function // here? +Matrix::Matrix (const boolMatrix& a) + : MArray2<double> (a.rows (), a.cols ()) +{ + for (int i = 0; i < a.rows (); i++) + for (int j = 0; j < a.cols (); j++) + elem (i, j) = a.elem (i, j); +} + Matrix::Matrix (const charMatrix& a) : MArray2<double> (a.rows (), a.cols ()) { @@ -1554,290 +1564,6 @@ return retval; } -// diagonal matrix by scalar -> matrix operations - -Matrix -operator + (const DiagMatrix& a, double s) -{ - Matrix tmp (a.rows (), a.cols (), s); - return a + tmp; -} - -Matrix -operator - (const DiagMatrix& a, double s) -{ - Matrix tmp (a.rows (), a.cols (), -s); - return a + tmp; -} - -// scalar by diagonal matrix -> matrix operations - -Matrix -operator + (double s, const DiagMatrix& a) -{ - Matrix tmp (a.rows (), a.cols (), s); - return tmp + a; -} - -Matrix -operator - (double s, const DiagMatrix& a) -{ - Matrix tmp (a.rows (), a.cols (), s); - return tmp - a; -} - -// matrix by diagonal matrix -> matrix operations - -Matrix -operator + (const Matrix& m, const DiagMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator +", nr, nc, a_nr, a_nc); - return Matrix (); - } - - if (nr == 0 || nc == 0) - return Matrix (nr, nc); - - Matrix result (m); - int a_len = a.length (); - for (int i = 0; i < a_len; i++) - result.elem (i, i) += a.elem (i, i); - - return result; -} - -Matrix -operator - (const Matrix& m, const DiagMatrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator -", nr, nc, a_nr, a_nc); - return Matrix (); - } - - if (nr == 0 || nc == 0) - return Matrix (nr, nc); - - Matrix result (m); - int a_len = a.length (); - for (int i = 0; i < a_len; i++) - result.elem (i, i) -= a.elem (i, i); - - return result; -} - -Matrix -operator * (const Matrix& m, const DiagMatrix& a) -{ - Matrix retval; - - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nc != a_nr) - gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); - else - { - if (nr == 0 || nc == 0 || a_nc == 0) - retval.resize (nr, a_nc, 0.0); - else - { - retval.resize (nr, a_nc); - double *c = retval.fortran_vec (); - - double *ctmp = 0; - - int a_len = a.length (); - - for (int j = 0; j < a_len; j++) - { - int idx = j * nr; - ctmp = c + idx; - - if (a.elem (j, j) == 1.0) - { - for (int i = 0; i < nr; i++) - ctmp[i] = m.elem (i, j); - } - else if (a.elem (j, j) == 0.0) - { - for (int i = 0; i < nr; i++) - ctmp[i] = 0.0; - } - else - { - for (int i = 0; i < nr; i++) - ctmp[i] = a.elem (j, j) * m.elem (i, j); - } - } - - if (a_nr < a_nc) - { - for (int i = nr * nc; i < nr * a_nc; i++) - ctmp[i] = 0.0; - } - } - } - - return retval; -} - -// diagonal matrix by matrix -> matrix operations - -Matrix -operator + (const DiagMatrix& m, const Matrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator +", nr, nc, a_nr, a_nc); - return Matrix (); - } - - if (nr == 0 || nc == 0) - return Matrix (nr, nc); - - Matrix result (a); - for (int i = 0; i < m.length (); i++) - result.elem (i, i) += m.elem (i, i); - - return result; -} - -Matrix -operator - (const DiagMatrix& m, const Matrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nr != a_nr || nc != a_nc) - { - gripe_nonconformant ("operator -", nr, nc, a_nr, a_nc); - return Matrix (); - } - - if (nr == 0 || nc == 0) - return Matrix (nr, nc); - - Matrix result (-a); - for (int i = 0; i < m.length (); i++) - result.elem (i, i) += m.elem (i, i); - - return result; -} - -Matrix -operator * (const DiagMatrix& m, const Matrix& a) -{ - int nr = m.rows (); - int nc = m.cols (); - int a_nr = a.rows (); - int a_nc = a.cols (); - if (nc != a_nr) - { - gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); - return Matrix (); - } - - if (nr == 0 || nc == 0 || a_nc == 0) - return Matrix (nr, a_nc, 0.0); - - Matrix c (nr, a_nc); - - for (int i = 0; i < m.length (); i++) - { - if (m.elem (i, i) == 1.0) - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = a.elem (i, j); - } - else if (m.elem (i, i) == 0.0) - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = 0.0; - } - else - { - for (int j = 0; j < a_nc; j++) - c.elem (i, j) = m.elem (i, i) * a.elem (i, j); - } - } - - if (nr > nc) - { - for (int j = 0; j < a_nc; j++) - for (int i = a_nr; i < nr; i++) - c.elem (i, j) = 0.0; - } - - return c; -} - -// matrix by matrix -> matrix operations - -Matrix -operator * (const Matrix& m, const Matrix& a) -{ - Matrix retval; - - int nr = m.rows (); - int nc = m.cols (); - - int a_nr = a.rows (); - int a_nc = a.cols (); - - if (nc != a_nr) - gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); - else - { - if (nr == 0 || nc == 0 || a_nc == 0) - retval.resize (nr, a_nc, 0.0); - else - { - int ld = nr; - int lda = a_nr; - - retval.resize (nr, a_nc); - double *c = retval.fortran_vec (); - - F77_XFCN (dgemm, DGEMM, ("N", "N", nr, a_nc, nc, 1.0, - m.data (), ld, a.data (), lda, 0.0, - c, nr, 1L, 1L)); - - if (f77_exception_encountered) - (*current_liboctave_error_handler) - ("unrecoverable error in dgemm"); - } - } - - return retval; -} - // other operations. Matrix @@ -3283,6 +3009,46 @@ return retval; } +// matrix by matrix -> matrix operations + +Matrix +operator * (const Matrix& m, const Matrix& a) +{ + Matrix retval; + + int nr = m.rows (); + int nc = m.cols (); + + int a_nr = a.rows (); + int a_nc = a.cols (); + + if (nc != a_nr) + gripe_nonconformant ("operator *", nr, nc, a_nr, a_nc); + else + { + if (nr == 0 || nc == 0 || a_nc == 0) + retval.resize (nr, a_nc, 0.0); + else + { + int ld = nr; + int lda = a_nr; + + retval.resize (nr, a_nc); + double *c = retval.fortran_vec (); + + F77_XFCN (dgemm, DGEMM, ("N", "N", nr, a_nc, nc, 1.0, + m.data (), ld, a.data (), lda, 0.0, + c, nr, 1L, 1L)); + + if (f77_exception_encountered) + (*current_liboctave_error_handler) + ("unrecoverable error in dgemm"); + } + } + + return retval; +} + /* ;;; Local Variables: *** ;;; mode: C++ ***
--- a/liboctave/dMatrix.h +++ b/liboctave/dMatrix.h @@ -60,6 +60,7 @@ // Matrix (const MDiagArray2<double>& a) : MArray2<double> (a) { } Matrix (const DiagMatrix& a); + Matrix (const boolMatrix& a); Matrix (const charMatrix& a); Matrix& operator = (const Matrix& a) @@ -175,32 +176,6 @@ friend Matrix operator * (const ColumnVector& a, const RowVector& b); - // diagonal matrix by scalar -> matrix operations - - friend Matrix operator + (const DiagMatrix& a, double s); - friend Matrix operator - (const DiagMatrix& a, double s); - - // scalar by diagonal matrix -> matrix operations - - friend Matrix operator + (double s, const DiagMatrix& a); - friend Matrix operator - (double s, const DiagMatrix& a); - - // matrix by diagonal matrix -> matrix operations - - friend Matrix operator + (const Matrix& a, const DiagMatrix& b); - friend Matrix operator - (const Matrix& a, const DiagMatrix& b); - friend Matrix operator * (const Matrix& a, const DiagMatrix& b); - - // diagonal matrix by matrix -> matrix operations - - friend Matrix operator + (const DiagMatrix& a, const Matrix& b); - friend Matrix operator - (const DiagMatrix& a, const Matrix& b); - friend Matrix operator * (const DiagMatrix& a, const Matrix& b); - - // matrix by matrix -> matrix operations - - friend Matrix operator * (const Matrix& a, const Matrix& b); - // other operations Matrix map (d_d_Mapper f) const; @@ -260,6 +235,8 @@ extern ComplexColumnVector Qzval (const Matrix& a, const Matrix& b); +extern Matrix operator * (const Matrix& a, const Matrix& b); + #endif /*
--- a/liboctave/idx-vector.cc +++ b/liboctave/idx-vector.cc @@ -33,6 +33,7 @@ #include <iostream.h> #include "Range.h" +#include "boolMatrix.h" #include "dColVector.h" #include "dMatrix.h" @@ -104,6 +105,7 @@ colon_equiv_checked = 0; colon_equiv = 0; colon = 0; + one_zero = 0; len = v.length (); @@ -114,7 +116,6 @@ { num_zeros = 0; num_ones = 0; - one_zero = 0; max_val = 0; min_val = 0; initialized = 1; @@ -146,6 +147,7 @@ colon_equiv_checked = 0; colon_equiv = 0; colon = 0; + one_zero = 0; orig_nr = m.rows (); orig_nc = m.columns (); @@ -156,7 +158,6 @@ { num_zeros = 0; num_ones = 0; - one_zero = 0; max_val = 0; min_val = 0; initialized = 1; @@ -190,6 +191,7 @@ colon_equiv_checked = 0; colon_equiv = 0; colon = 0; + one_zero = 0; len = 1; @@ -216,6 +218,7 @@ colon_equiv_checked = 0; colon_equiv = 0; colon = 0; + one_zero = 0; len = r.nelem (); @@ -231,7 +234,6 @@ { num_zeros = 0; num_ones = 0; - one_zero = 0; max_val = 0; min_val = 0; initialized = 1; @@ -274,6 +276,66 @@ init_state (); } +IDX_VEC_REP::idx_vector_rep (bool b) +{ + data = 0; + initialized = 0; + frozen = 0; + colon_equiv_checked = 0; + colon_equiv = 0; + colon = 0; + one_zero = 1; + + len = 1; + + orig_nr = 1; + orig_nc = 1; + + data = new int [len]; + + data[0] = tree_to_mat_idx (b); + + init_state (); +} + +IDX_VEC_REP::idx_vector_rep (const boolMatrix& bm) +{ + data = 0; + initialized = 0; + frozen = 0; + colon_equiv_checked = 0; + colon_equiv = 0; + colon = 0; + one_zero = 1; + + orig_nr = bm.rows (); + orig_nc = bm.columns (); + + len = orig_nr * orig_nc; + + if (len == 0) + { + num_zeros = 0; + num_ones = 0; + one_zero = 0; + max_val = 0; + min_val = 0; + initialized = 1; + return; + } + else + { + int k = 0; + data = new int [len]; + + for (int j = 0; j < orig_nc; j++) + for (int i = 0; i < orig_nr; i++) + data[k++] = tree_to_mat_idx (bm.elem (i, j)); + } + + init_state (); +} + IDX_VEC_REP& IDX_VEC_REP::operator = (const IDX_VEC_REP& a) { @@ -313,13 +375,11 @@ if (colon) { - one_zero = 0; - min_val = max_val = 0; + min_val = 0; + max_val = 0; } else { - one_zero = 1; - min_val = max_val = data[0]; int i = 0; @@ -330,9 +390,6 @@ else if (data[i] == 0) num_ones++; - if (one_zero && data[i] != -1 && data[i] != 0) - one_zero = 0; - if (data[i] > max_val) max_val = data[i]; @@ -346,10 +403,9 @@ } void -IDX_VEC_REP::maybe_convert_one_zero_to_idx (int z_len, int prefer_zero_one) +IDX_VEC_REP::maybe_convert_one_zero_to_idx (int z_len) { - if (one_zero && z_len == len - && (num_ones != len || prefer_zero_one)) + if (one_zero && z_len == len) { if (num_ones == 0) { @@ -515,7 +571,7 @@ frozen_len = 0; else { - maybe_convert_one_zero_to_idx (z_len, prefer_zero_one); + maybe_convert_one_zero_to_idx (z_len); max_val = max (); min_val = min ();
--- a/liboctave/idx-vector.h +++ b/liboctave/idx-vector.h @@ -29,6 +29,7 @@ class ostream; class ColumnVector; +class boolMatrix; class Matrix; class Range; @@ -68,6 +69,10 @@ idx_vector_rep (char c); + idx_vector_rep (bool b); + + idx_vector_rep (const boolMatrix& bm); + idx_vector_rep (const idx_vector_rep& a); ~idx_vector_rep (void) { delete [] data; } @@ -128,7 +133,7 @@ void init_state (void); - void maybe_convert_one_zero_to_idx (int z_len, int prefer_zero_one); + void maybe_convert_one_zero_to_idx (int z_len); }; public: @@ -169,6 +174,18 @@ rep->count = 1; } + idx_vector (bool b) + { + rep = new idx_vector_rep (b); + rep->count = 1; + } + + idx_vector (const boolMatrix& bm) + { + rep = new idx_vector_rep (bm); + rep->count = 1; + } + idx_vector (const idx_vector& a) { rep = a.rep; @@ -238,8 +255,8 @@ friend ostream& operator << (ostream& os, const idx_vector& a) { return a.print (os); } - void maybe_convert_one_zero_to_idx (int z_len, int prefer_zero_one = 0) - { rep->maybe_convert_one_zero_to_idx (z_len, prefer_zero_one); } + void maybe_convert_one_zero_to_idx (int z_len) + { rep->maybe_convert_one_zero_to_idx (z_len); } private:
--- a/liboctave/mx-base.h +++ b/liboctave/mx-base.h @@ -25,6 +25,7 @@ // Matrix classes. +#include "boolMatrix.h" #include "chMatrix.h" #include "dMatrix.h" #include "CMatrix.h"
--- a/liboctave/mx-defs.h +++ b/liboctave/mx-defs.h @@ -27,6 +27,7 @@ class Matrix; class ComplexMatrix; +class boolMatrix; class charMatrix; class ColumnVector;
--- a/liboctave/mx-inlines.cc +++ b/liboctave/mx-inlines.cc @@ -20,8 +20,8 @@ */ -#if !defined (octave_mx_ops_h) -#define octave_mx_ops 1 +#if !defined (octave_mx_inlines_h) +#define octave_mx_inlines_h 1 #include <cstddef> @@ -194,6 +194,7 @@ return true; \ } +OP_EQ_FCN (bool, bool) OP_EQ_FCN (char, char) OP_EQ_FCN (double, double) OP_EQ_FCN (Complex, Complex)