Mercurial > hg > octave-lojdl
changeset 8375:e3c9102431a9
fix design problems of diag & perm matrix classes
author | Jaroslav Hajek <highegg@gmail.com> |
---|---|
date | Fri, 05 Dec 2008 10:20:18 +0100 |
parents | 0f8e810c6950 |
children | c43481a19bfe |
files | liboctave/CDiagMatrix.cc liboctave/CDiagMatrix.h liboctave/ChangeLog liboctave/CmplxQRP.cc liboctave/DiagArray2.cc liboctave/DiagArray2.h liboctave/MDiagArray2.h liboctave/PermMatrix.cc liboctave/PermMatrix.h liboctave/dDiagMatrix.cc liboctave/dDiagMatrix.h liboctave/dMatrix.cc liboctave/dbleQRP.cc liboctave/fCDiagMatrix.cc liboctave/fCDiagMatrix.h liboctave/fCmplxQRP.cc liboctave/fDiagMatrix.cc liboctave/fDiagMatrix.h liboctave/fMatrix.cc liboctave/floatQRP.cc liboctave/idx-vector.cc liboctave/mx-op-defs.h |
diffstat | 22 files changed, 155 insertions(+), 204 deletions(-) [+] |
line wrap: on
line diff
--- a/liboctave/CDiagMatrix.cc +++ b/liboctave/CDiagMatrix.cc @@ -529,47 +529,6 @@ // other operations -ComplexColumnVector -ComplexDiagMatrix::diag (octave_idx_type k) const -{ - octave_idx_type nnr = rows (); - octave_idx_type nnc = cols (); - if (k > 0) - nnc -= k; - else if (k < 0) - nnr += k; - - ComplexColumnVector d; - - if (nnr > 0 && nnc > 0) - { - octave_idx_type ndiag = (nnr < nnc) ? nnr : nnc; - - d.resize (ndiag); - - if (k > 0) - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i, i+k); - } - else if ( k < 0) - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i-k, i); - } - else - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i, i); - } - } - else - (*current_liboctave_error_handler) - ("diag: requested diagonal out of range"); - - return d; -} - ComplexDET ComplexDiagMatrix::determinant (void) const {
--- a/liboctave/CDiagMatrix.h +++ b/liboctave/CDiagMatrix.h @@ -122,7 +122,8 @@ // other operations - ComplexColumnVector diag (octave_idx_type k = 0) const; + ComplexColumnVector diag (octave_idx_type k = 0) const + { return MDiagArray2<Complex>::diag (k); } ComplexDET determinant (void) const; double rcond (void) const;
--- a/liboctave/ChangeLog +++ b/liboctave/ChangeLog @@ -1,3 +1,21 @@ +2008-12-04 Jaroslav Hajek <highegg@gmail.com> + + * DiagArray2.h (DiagArray2<T>): Inherit Array<T> privately. + (DiagArray2<T>::dim1, dim2, rows, columns, cols, length, + nelem, numel, byte_size, dims): New methods. + (DiagArray2<T>::diag): New method decl. + * DiagArray2.cc (DiagArray2<T>::diag): New method. + * MDiagArray2.h (MDiagArray2<T>::diag): New method. + * dDiagMatrix.cc (DiagMatrix::diag): Remove. + * fDiagMatrix.cc (FloatDiagMatrix::diag): Remove. + * CDiagMatrix.cc (ComplexDiagMatrix::diag): Remove. + * fCDiagMatrix.cc (FloatComplexDiagMatrix::diag): Remove. + + * PermMatrix.h (PermMatrix): Inherit Array<octave_idx_type> privately. + (PermMatrix::dim1, dim2, rows, columns, cols, length, + nelem, numel, byte_size, dims): New methods. + + 2008-12-04 Jaroslav Hajek <highegg@gmail.com> * dDiagMatrix.cc (DiagMatrix::determinant, DiagMatrix::rcond): New
--- a/liboctave/CmplxQRP.cc +++ b/liboctave/CmplxQRP.cc @@ -129,7 +129,7 @@ ColumnVector ComplexQRP::Pvec (void) const { - Array<double> pa (p); + Array<double> pa (p.pvec ()); ColumnVector pv (MArray<double> (pa) + 1.0); return pv; }
--- a/liboctave/DiagArray2.cc +++ b/liboctave/DiagArray2.cc @@ -30,16 +30,41 @@ #include <iostream> +#include <algorithm> + #include "DiagArray2.h" #include "lo-error.h" template <class T> +Array<T> +DiagArray2<T>::diag (octave_idx_type k) const +{ + Array<T> d; + + if (k == 0) + { + // The main diagonal is shallow-copied. + d = *this; + d.dimensions = dim_vector (length ()); + } + else if (k > 0 && k < cols ()) + d = Array<T> (std::min (cols () - k, rows ()), T ()); + else if (k < 0 && -k < rows ()) + d = Array<T> (std::min (rows () + k, cols ()), T ()); + else + (*current_liboctave_error_handler) + ("diag: requested diagonal out of range"); + + return d; +} + +template <class T> DiagArray2<T> DiagArray2<T>::transpose (void) const { DiagArray2<T> retval (*this); - retval.dimensions = dim_vector (this->dim2 (), this->dim1 ()); + retval.dimensions = dim_vector (dim2 (), dim1 ()); return retval; } @@ -47,7 +72,7 @@ DiagArray2<T> DiagArray2<T>::hermitian (T (* fcn) (const T&)) const { - DiagArray2<T> retval (this->dim2 (), this->dim1 ()); + DiagArray2<T> retval (dim2 (), dim1 ()); const T *p = this->data (); T *q = retval.fortran_vec (); for (octave_idx_type i = 0; i < this->length (); i++) @@ -61,7 +86,7 @@ T DiagArray2<T>::checkelem (octave_idx_type r, octave_idx_type c) const { - if (r < 0 || c < 0 || r >= this->dim1 () || c >= this->dim2 ()) + if (r < 0 || c < 0 || r >= dim1 () || c >= dim2 ()) { (*current_liboctave_error_handler) ("range error in DiagArray2"); return T (); @@ -79,7 +104,7 @@ return; } - if (r == this->dim1 () && c == this->dim2 ()) + if (r == dim1 () && c == dim2 ()) return; typename Array<T>::ArrayRep *old_rep = Array<T>::rep; @@ -114,7 +139,7 @@ return; } - if (r == this->dim1 () && c == this->dim2 ()) + if (r == dim1 () && c == dim2 ()) return; typename Array<T>::ArrayRep *old_rep = Array<T>::rep;
--- a/liboctave/DiagArray2.h +++ b/liboctave/DiagArray2.h @@ -3,6 +3,7 @@ Copyright (C) 1996, 1997, 2000, 2002, 2003, 2004, 2005, 2006, 2007 John W. Eaton +Copyright (C) 2008 Jaroslav Hajek This file is part of Octave. @@ -31,8 +32,6 @@ #include "Array.h" #include "lo-error.h" -class idx_vector; - // A two-dimensional array with diagonal elements only. // // Idea and example code for Proxy class and functions from: @@ -46,9 +45,13 @@ // James Kanze email: kanze@us-es.sel.de // GABI Software, Sarl., 8 rue du Faisan, F-67000 Strasbourg, France +// Array<T> is inherited privately because we abuse the dimensions variable +// for true dimensions. Therefore, the inherited Array<T> object is not a valid +// Array<T> object, and should not be publicly accessible. + template <class T> class -DiagArray2 : public Array<T> +DiagArray2 : private Array<T> { private: @@ -112,6 +115,8 @@ public: + typedef T element_type; + DiagArray2 (void) : Array<T> (dim_vector (0, 0)) { } DiagArray2 (octave_idx_type r, octave_idx_type c) : Array<T> (r < c ? r : c) @@ -131,7 +136,7 @@ { this->dimensions = a.dims (); } template <class U> - DiagArray2 (const DiagArray2<U>& a) : Array<T> (a) + DiagArray2 (const DiagArray2<U>& a) : Array<T> (a.diag ()) { this->dimensions = a.dims (); } ~DiagArray2 (void) { } @@ -144,6 +149,24 @@ return *this; } + + octave_idx_type dim1 (void) const { return Array<T>::dimensions(0); } + octave_idx_type dim2 (void) const { return Array<T>::dimensions(1); } + + octave_idx_type rows (void) const { return dim1 (); } + octave_idx_type cols (void) const { return dim2 (); } + octave_idx_type columns (void) const { return dim2 (); } + + octave_idx_type length (void) const { return Array<T>::length (); } + octave_idx_type nelem (void) const { return dim1 () * dim2 (); } + octave_idx_type numel (void) const { return nelem (); } + + size_t byte_size (void) const { return length () * sizeof (T); } + + dim_vector dims (void) const { return Array<T>::dimensions; } + + Array<T> diag (octave_idx_type k = 0) const; + Proxy elem (octave_idx_type r, octave_idx_type c) { return Proxy (this, r, c); @@ -151,7 +174,7 @@ Proxy checkelem (octave_idx_type r, octave_idx_type c) { - if (r < 0 || c < 0 || r >= this->dim1 () || c >= this->dim2 ()) + if (r < 0 || c < 0 || r >= dim1 () || c >= dim2 ()) { (*current_liboctave_error_handler) ("range error in DiagArray2"); return Proxy (0, r, c); @@ -162,7 +185,7 @@ Proxy operator () (octave_idx_type r, octave_idx_type c) { - if (r < 0 || c < 0 || r >= this->dim1 () || c >= this->dim2 ()) + if (r < 0 || c < 0 || r >= dim1 () || c >= dim2 ()) { (*current_liboctave_error_handler) ("range error in DiagArray2"); return Proxy (0, r, c); @@ -204,6 +227,15 @@ DiagArray2<T> transpose (void) const; DiagArray2<T> hermitian (T (*fcn) (const T&) = 0) const; + + const T *data (void) const { return Array<T>::data (); } + + const T *fortran_vec (void) const { return Array<T>::fortran_vec (); } + + T *fortran_vec (void) { return Array<T>::fortran_vec (); } + + void print_info (std::ostream& os, const std::string& prefix) const + { Array<T>::print_info (os, prefix); } }; #endif
--- a/liboctave/MDiagArray2.h +++ b/liboctave/MDiagArray2.h @@ -27,6 +27,7 @@ #include "DiagArray2.h" #include "MArray2.h" +#include "MArray.h" // Two dimensional diagonal array with math ops. @@ -88,9 +89,9 @@ { octave_idx_type retval = 0; - const T *d = this->Array<T>::data (); + const T *d = this->data (); - octave_idx_type nel = this->Array<T>::numel (); + octave_idx_type nel = this->length (); for (octave_idx_type i = 0; i < nel; i++) { @@ -101,6 +102,9 @@ return retval; } + MArray<T> diag (octave_idx_type k = 0) const + { return DiagArray2<T>::diag (k); } + MDiagArray2<T> transpose (void) const { return DiagArray2<T>::transpose (); } MDiagArray2<T> hermitian (T (*fcn) (const T&) = 0) const { return DiagArray2<T>::hermitian (fcn); }
--- a/liboctave/PermMatrix.cc +++ b/liboctave/PermMatrix.cc @@ -123,7 +123,7 @@ PermMatrix operator *(const PermMatrix& a, const PermMatrix& b) { - const Array<octave_idx_type>& ia = a, ib = b; + const Array<octave_idx_type> ia = a.pvec (), ib = b.pvec (); PermMatrix r; octave_idx_type n = a.columns (); if (n != b.rows ())
--- a/liboctave/PermMatrix.h +++ b/liboctave/PermMatrix.h @@ -26,7 +26,11 @@ #include "Array.h" #include "mx-defs.h" -class PermMatrix : public Array<octave_idx_type> +// Array<T> is inherited privately because we abuse the dimensions variable +// for true dimensions. Therefore, the inherited Array<T> object is not a valid +// Array<T> object, and should not be publicly accessible. + +class PermMatrix : private Array<octave_idx_type> { private: @@ -42,11 +46,35 @@ bool check = true); PermMatrix (const PermMatrix& m) - : Array<octave_idx_type> (m), _colp(m._colp) - { this->dimensions = m.dims (); } + : Array<octave_idx_type> (m), _colp(m._colp) { } PermMatrix (const idx_vector& idx, bool colp = false, octave_idx_type n = 0); + octave_idx_type dim1 (void) const + { return Array<octave_idx_type>::dimensions(0); } + octave_idx_type dim2 (void) const + { return Array<octave_idx_type>::dimensions(1); } + + octave_idx_type rows (void) const { return dim1 (); } + octave_idx_type cols (void) const { return dim2 (); } + octave_idx_type columns (void) const { return dim2 (); } + + octave_idx_type length (void) const + { return Array<octave_idx_type>::length (); } + octave_idx_type nelem (void) const { return dim1 () * dim2 (); } + octave_idx_type numel (void) const { return nelem (); } + + size_t byte_size (void) const { return length () * sizeof (octave_idx_type); } + + dim_vector dims (void) const { return Array<octave_idx_type>::dimensions; } + + Array<octave_idx_type> pvec (void) const + { + Array<octave_idx_type> retval (*this); + retval.dimensions = dim_vector (length ()); + return retval; + } + octave_idx_type elem (octave_idx_type i, octave_idx_type j) const { @@ -80,6 +108,18 @@ friend PermMatrix operator *(const PermMatrix& a, const PermMatrix& b); + const octave_idx_type *data (void) const + { return Array<octave_idx_type>::data (); } + + const octave_idx_type *fortran_vec (void) const + { return Array<octave_idx_type>::fortran_vec (); } + + octave_idx_type *fortran_vec (void) + { return Array<octave_idx_type>::fortran_vec (); } + + void print_info (std::ostream& os, const std::string& prefix) const + { Array<octave_idx_type>::print_info (os, prefix); } + private: bool _colp; };
--- a/liboctave/dDiagMatrix.cc +++ b/liboctave/dDiagMatrix.cc @@ -342,51 +342,6 @@ // other operations -ColumnVector -DiagMatrix::diag (octave_idx_type k) const -{ - ColumnVector d; - - octave_idx_type nnr = rows (); - octave_idx_type nnc = cols (); - - if (nnr == 0 || nnc == 0) - return d; - - if (k > 0) - nnc -= k; - else if (k < 0) - nnr += k; - - if (nnr > 0 && nnc > 0) - { - octave_idx_type ndiag = (nnr < nnc) ? nnr : nnc; - - d.resize (ndiag); - - if (k > 0) - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i, i+k); - } - else if ( k < 0) - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i-k, i); - } - else - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i, i); - } - } - else - (*current_liboctave_error_handler) - ("diag: requested diagonal out of range"); - - return d; -} - DET DiagMatrix::determinant (void) const {
--- a/liboctave/dDiagMatrix.h +++ b/liboctave/dDiagMatrix.h @@ -97,7 +97,8 @@ // other operations - ColumnVector diag (octave_idx_type k = 0) const; + ColumnVector diag (octave_idx_type k = 0) const + { return MDiagArray2<double>::diag (k); } DET determinant (void) const; double rcond (void) const;
--- a/liboctave/dMatrix.cc +++ b/liboctave/dMatrix.cc @@ -251,7 +251,7 @@ Matrix::Matrix (const PermMatrix& a) : MArray2<double> (a.rows (), a.cols (), 0.0) { - const Array<octave_idx_type> ia (a); + const Array<octave_idx_type> ia (a.pvec ()); octave_idx_type len = a.rows (); if (a.is_col_perm ()) for (octave_idx_type i = 0; i < len; i++)
--- a/liboctave/dbleQRP.cc +++ b/liboctave/dbleQRP.cc @@ -123,7 +123,7 @@ ColumnVector QRP::Pvec (void) const { - Array<double> pa (p); + Array<double> pa (p.pvec ()); ColumnVector pv (MArray<double> (pa) + 1.0); return pv; }
--- a/liboctave/fCDiagMatrix.cc +++ b/liboctave/fCDiagMatrix.cc @@ -550,47 +550,6 @@ // other operations -FloatComplexColumnVector -FloatComplexDiagMatrix::diag (octave_idx_type k) const -{ - octave_idx_type nnr = rows (); - octave_idx_type nnc = cols (); - if (k > 0) - nnc -= k; - else if (k < 0) - nnr += k; - - FloatComplexColumnVector d; - - if (nnr > 0 && nnc > 0) - { - octave_idx_type ndiag = (nnr < nnc) ? nnr : nnc; - - d.resize (ndiag); - - if (k > 0) - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i, i+k); - } - else if ( k < 0) - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i-k, i); - } - else - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i, i); - } - } - else - (*current_liboctave_error_handler) - ("diag: requested diagonal out of range"); - - return d; -} - FloatComplexDET FloatComplexDiagMatrix::determinant (void) const {
--- a/liboctave/fCDiagMatrix.h +++ b/liboctave/fCDiagMatrix.h @@ -122,7 +122,8 @@ // other operations - FloatComplexColumnVector diag (octave_idx_type k = 0) const; + FloatComplexColumnVector diag (octave_idx_type k = 0) const + { return MDiagArray2<FloatComplex>::diag (k); } FloatComplexDET determinant (void) const; float rcond (void) const;
--- a/liboctave/fCmplxQRP.cc +++ b/liboctave/fCmplxQRP.cc @@ -129,7 +129,7 @@ FloatColumnVector FloatComplexQRP::Pvec (void) const { - Array<float> pa (p); + Array<float> pa (p.pvec ()); FloatColumnVector pv (MArray<float> (pa) + 1.0f); return pv; }
--- a/liboctave/fDiagMatrix.cc +++ b/liboctave/fDiagMatrix.cc @@ -349,51 +349,6 @@ // other operations -FloatColumnVector -FloatDiagMatrix::diag (octave_idx_type k) const -{ - FloatColumnVector d; - - octave_idx_type nnr = rows (); - octave_idx_type nnc = cols (); - - if (nnr == 0 || nnc == 0) - return d; - - if (k > 0) - nnc -= k; - else if (k < 0) - nnr += k; - - if (nnr > 0 && nnc > 0) - { - octave_idx_type ndiag = (nnr < nnc) ? nnr : nnc; - - d.resize (ndiag); - - if (k > 0) - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i, i+k); - } - else if ( k < 0) - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i-k, i); - } - else - { - for (octave_idx_type i = 0; i < ndiag; i++) - d.elem (i) = elem (i, i); - } - } - else - (*current_liboctave_error_handler) - ("diag: requested diagonal out of range"); - - return d; -} - FloatDET FloatDiagMatrix::determinant (void) const {
--- a/liboctave/fDiagMatrix.h +++ b/liboctave/fDiagMatrix.h @@ -97,7 +97,8 @@ // other operations - FloatColumnVector diag (octave_idx_type k = 0) const; + FloatColumnVector diag (octave_idx_type k = 0) const + { return MDiagArray2<float>::diag (k); } FloatDET determinant (void) const; float rcond (void) const;
--- a/liboctave/fMatrix.cc +++ b/liboctave/fMatrix.cc @@ -250,7 +250,7 @@ FloatMatrix::FloatMatrix (const PermMatrix& a) : MArray2<float> (a.rows (), a.cols (), 0.0) { - const Array<octave_idx_type> ia (a); + const Array<octave_idx_type> ia (a.pvec ()); octave_idx_type len = a.rows (); if (a.is_col_perm ()) for (octave_idx_type i = 0; i < len; i++)
--- a/liboctave/floatQRP.cc +++ b/liboctave/floatQRP.cc @@ -123,7 +123,7 @@ FloatColumnVector FloatQRP::Pvec (void) const { - Array<float> pa (p); + Array<float> pa (p.pvec ()); FloatColumnVector pv (MArray<float> (pa) + 1.0f); return pv; }
--- a/liboctave/idx-vector.cc +++ b/liboctave/idx-vector.cc @@ -557,7 +557,7 @@ } } - delete left; + delete [] left; } return retval;
--- a/liboctave/mx-op-defs.h +++ b/liboctave/mx-op-defs.h @@ -1301,10 +1301,10 @@ if (p.is_col_perm ()) \ { \ result = M (nr, nc); \ - result.assign (idx_vector (p), idx_vector::colon, x); \ + result.assign (p.pvec (), idx_vector::colon, x); \ } \ else \ - result = x.index (idx_vector (p), idx_vector::colon); \ + result = x.index (p.pvec (), idx_vector::colon); \ } \ \ return result; \ @@ -1320,11 +1320,11 @@ else \ { \ if (p.is_col_perm ()) \ - result = x.index (idx_vector::colon, idx_vector (p)); \ + result = x.index (idx_vector::colon, p.pvec ()); \ else \ { \ result = M (nr, nc); \ - result.assign (idx_vector::colon, idx_vector (p), x); \ + result.assign (idx_vector::colon, p.pvec (), x); \ } \ } \ \