Mercurial > hg > octave-terminal
changeset 4725:fa612b2cbfe9
[project @ 2004-01-23 16:42:51 by jwe]
author | jwe |
---|---|
date | Fri, 23 Jan 2004 16:42:51 +0000 |
parents | bdacd0383fbd |
children | 14dc2267c343 |
files | liboctave/ChangeLog liboctave/EIG.cc liboctave/EIG.h src/ChangeLog src/DLD-FUNCTIONS/eig.cc |
diffstat | 5 files changed, 46 insertions(+), 27 deletions(-) [+] |
line wrap: on
line diff
--- a/liboctave/ChangeLog +++ b/liboctave/ChangeLog @@ -1,3 +1,10 @@ +2004-01-23 John W. Eaton <jwe@bevo.che.wisc.edu> + + * EIG.cc, EIG.h (EIG::init, EIG::symmetric_init, EIG::hermitian_init): + New arg, calc_eigenvectors. + * EIG.h (EIG:EIG): New optional arg, calc_eigenvectors. + Based on patch from David Bateman <dbateman@free.fr>. + 2004-01-22 John W. Eaton <jwe@bevo.che.wisc.edu> * Array.cc (Array<T>::assign2, Array<T>::assignN):
--- a/liboctave/EIG.cc +++ b/liboctave/EIG.cc @@ -71,10 +71,10 @@ } int -EIG::init (const Matrix& a) +EIG::init (const Matrix& a, bool calc_ev) { if (a.is_symmetric ()) - return symmetric_init (a); + return symmetric_init (a, calc_ev); int n = a.rows (); @@ -95,7 +95,8 @@ Array<double> wi (n); double *pwi = wi.fortran_vec (); - Matrix vr (n, n); + int nvr = calc_ev ? n : 0; + Matrix vr (nvr, nvr); double *pvr = vr.fortran_vec (); // XXX FIXME XXX -- it might be possible to choose a better value of @@ -109,7 +110,7 @@ int idummy = 1; F77_XFCN (dgeev, DGEEV, (F77_CONST_CHAR_ARG2 ("N", 1), - F77_CONST_CHAR_ARG2 ("V", 1), + F77_CONST_CHAR_ARG2 (calc_ev ? "V" : "N", 1), n, tmp_data, n, pwr, pwi, dummy, idummy, pvr, n, pwork, lwork, info F77_CHAR_ARG_LEN (1) @@ -124,14 +125,14 @@ else { lambda.resize (n); - v.resize (n, n); + v.resize (nvr, nvr); for (int j = 0; j < n; j++) { if (wi.elem (j) == 0.0) { lambda.elem (j) = Complex (wr.elem (j)); - for (int i = 0; i < n; i++) + for (int i = 0; i < nvr; i++) v.elem (i, j) = vr.elem (i, j); } else @@ -146,7 +147,7 @@ lambda.elem(j) = Complex (wr.elem(j), wi.elem(j)); lambda.elem(j+1) = Complex (wr.elem(j+1), wi.elem(j+1)); - for (int i = 0; i < n; i++) + for (int i = 0; i < nvr; i++) { double real_part = vr.elem (i, j); double imag_part = vr.elem (i, j+1); @@ -163,7 +164,7 @@ } int -EIG::symmetric_init (const Matrix& a) +EIG::symmetric_init (const Matrix& a, bool calc_ev) { int n = a.rows (); @@ -188,7 +189,7 @@ Array<double> work (lwork); double *pwork = work.fortran_vec (); - F77_XFCN (dsyev, DSYEV, (F77_CONST_CHAR_ARG2 ("V", 1), + F77_XFCN (dsyev, DSYEV, (F77_CONST_CHAR_ARG2 (calc_ev ? "V" : "N", 1), F77_CONST_CHAR_ARG2 ("U", 1), n, tmp_data, n, pwr, pwork, lwork, info F77_CHAR_ARG_LEN (1) @@ -201,17 +202,17 @@ else { lambda = ComplexColumnVector (wr); - v = ComplexMatrix (atmp); + v = calc_ev ? ComplexMatrix (atmp) : ComplexMatrix (); } return info; } int -EIG::init (const ComplexMatrix& a) +EIG::init (const ComplexMatrix& a, bool calc_ev) { if (a.is_hermitian ()) - return hermitian_init (a); + return hermitian_init (a, calc_ev); int n = a.rows (); @@ -229,7 +230,8 @@ ComplexColumnVector w (n); Complex *pw = w.fortran_vec (); - ComplexMatrix vtmp (n, n); + int nvr = calc_ev ? n : 0; + ComplexMatrix vtmp (nvr, nvr); Complex *pv = vtmp.fortran_vec (); // XXX FIXME XXX -- it might be possible to choose a better value of @@ -247,7 +249,7 @@ int idummy = 1; F77_XFCN (zgeev, ZGEEV, (F77_CONST_CHAR_ARG2 ("N", 1), - F77_CONST_CHAR_ARG2 ("V", 1), + F77_CONST_CHAR_ARG2 (calc_ev ? "V" : "N", 1), n, tmp_data, n, pw, dummy, idummy, pv, n, pwork, lwork, prwork, info F77_CHAR_ARG_LEN (1) @@ -267,7 +269,7 @@ } int -EIG::hermitian_init (const ComplexMatrix& a) +EIG::hermitian_init (const ComplexMatrix& a, bool calc_ev) { int n = a.rows (); @@ -296,7 +298,7 @@ Array<double> rwork (lrwork); double *prwork = rwork.fortran_vec (); - F77_XFCN (zheev, ZHEEV, (F77_CONST_CHAR_ARG2 ("V", 1), + F77_XFCN (zheev, ZHEEV, (F77_CONST_CHAR_ARG2 (calc_ev ? "V" : "N", 1), F77_CONST_CHAR_ARG2 ("U", 1), n, tmp_data, n, pwr, pwork, lwork, prwork, info F77_CHAR_ARG_LEN (1) @@ -309,7 +311,7 @@ else { lambda = ComplexColumnVector (wr); - v = ComplexMatrix (atmp); + v = calc_ev ? ComplexMatrix (atmp) : ComplexMatrix (); } return info;
--- a/liboctave/EIG.h +++ b/liboctave/EIG.h @@ -44,13 +44,17 @@ EIG (void) : lambda (), v () { } - EIG (const Matrix& a) { init (a); } + EIG (const Matrix& a, bool calc_eigenvectors = true) + { init (a, calc_eigenvectors); } - EIG (const Matrix& a, int& info) { info = init (a); } + EIG (const Matrix& a, int& info, bool calc_eigenvectors = true) + { info = init (a, calc_eigenvectors); } - EIG (const ComplexMatrix& a) { init (a); } + EIG (const ComplexMatrix& a, bool calc_eigenvectors = true) + { init (a, calc_eigenvectors); } - EIG (const ComplexMatrix& a, int& info) { info = init (a); } + EIG (const ComplexMatrix& a, int& info, bool calc_eigenvectors = true) + { info = init (a, calc_eigenvectors); } EIG (const EIG& a) : lambda (a.lambda), v (a.v) { } @@ -78,11 +82,11 @@ ComplexColumnVector lambda; ComplexMatrix v; - int init (const Matrix& a); - int init (const ComplexMatrix& a); + int init (const Matrix& a, bool calc_eigenvectors); + int init (const ComplexMatrix& a, bool calc_eigenvectors); - int symmetric_init (const Matrix& a); - int hermitian_init (const ComplexMatrix& a); + int symmetric_init (const Matrix& a, bool calc_eigenvectors); + int hermitian_init (const ComplexMatrix& a, bool calc_eigenvectors); }; #endif
--- a/src/ChangeLog +++ b/src/ChangeLog @@ -1,3 +1,9 @@ +2004-01-23 John W. Eaton <jwe@bevo.che.wisc.edu> + + * DLD-FUNCTIONS/eig.cc (Feig): Use new optional arg for EIG to + avoid computing eigenvectors if not requested. + Based on a patch from David Bateman <dbateman@free.fr>. + 2004-01-23 John W. Eaton <jwe@bevo.che.wisc.edu> * ov-cell.cc (all_strings): Always compute total required length
--- a/src/DLD-FUNCTIONS/eig.cc +++ b/src/DLD-FUNCTIONS/eig.cc @@ -81,7 +81,7 @@ if (error_state) return retval; else - result = EIG (tmp); + result = EIG (tmp, nargout > 1); } else if (arg.is_complex_type ()) { @@ -90,7 +90,7 @@ if (error_state) return retval; else - result = EIG (ctmp); + result = EIG (ctmp, nargout > 1); } else {