octave: liboctave/numeric/svd.cc annotate

annotate liboctave/numeric/svd.cc @ 25349:a8a4b6e6e754 stable

Fix lrwork dimension in svd for complex cases (bug #53823). * svd.cc: Due to a bug in LAPACK <= 3.6, if jobz == 'N' then lrwork must be set to 7*min_mn.

author	Marco Caliari <marco.caliari@univr.it>
date	Thu, 03 May 2018 16:47:30 +0200
parents	6652d3823428
children	00f796120a6d

rev	line source
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	1 /*
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	2
25054 6652d3823428 maint: Update copyright dates in all source files. John W. Eaton <jwe@octave.org> parents: 24534 diff changeset	3 Copyright (C) 2016-2018 Carnë Draug
6652d3823428 maint: Update copyright dates in all source files. John W. Eaton <jwe@octave.org> parents: 24534 diff changeset	4 Copyright (C) 1994-2018 John W. Eaton
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	5
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	6 This file is part of Octave.
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	7
24534 194eb4bd202b maint: Update punctuation for GPL v3 license text. Rik <rik@octave.org> parents: 23696 diff changeset	8 Octave is free software: you can redistribute it and/or modify it
22755 3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22592 diff changeset	9 under the terms of the GNU General Public License as published by
24534 194eb4bd202b maint: Update punctuation for GPL v3 license text. Rik <rik@octave.org> parents: 23696 diff changeset	10 the Free Software Foundation, either version 3 of the License, or
22755 3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22592 diff changeset	11 (at your option) any later version.
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	12
22755 3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22592 diff changeset	13 Octave is distributed in the hope that it will be useful, but
3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22592 diff changeset	14 WITHOUT ANY WARRANTY; without even the implied warranty of
3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22592 diff changeset	15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
3a2b891d0b33 maint: Standardize Copyright formatting. Rik <rik@octave.org> parents: 22592 diff changeset	16 GNU General Public License for more details.
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	17
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	18 You should have received a copy of the GNU General Public License
7016 93c65f2a5668 [project @ 2007-10-12 06:40:56 by jwe] jwe parents: 5307 diff changeset	19 along with Octave; see the file COPYING. If not, see
24534 194eb4bd202b maint: Update punctuation for GPL v3 license text. Rik <rik@octave.org> parents: 23696 diff changeset	20 <https://www.gnu.org/licenses/>.
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	21
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	22 */
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	23
21724 aba2e6293dd8 use "#if ..." consistently instead of "#ifdef" and "#ifndef" John W. Eaton <jwe@octave.org> parents: 21301 diff changeset	24 #if defined (HAVE_CONFIG_H)
21301 40de9f8f23a6 Use '#include "config.h"' rather than <config.h>. Rik <rik@octave.org> parents: 21273 diff changeset	25 # include "config.h"
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	26 #endif
3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	27
23475 d691ed308237 maint: Clean up #includes in liboctave/numeric directory. Rik <rik@octave.org> parents: 23450 diff changeset	28 #include <cassert>
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	29
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	30 #include <algorithm>
1631 b0b22b6ce22f [project @ 1995-11-14 20:24:15 by jwe] jwe parents: 1545 diff changeset	31
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	32 #include "CMatrix.h"
cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	33 #include "dDiagMatrix.h"
cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	34 #include "dMatrix.h"
cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	35 #include "fCMatrix.h"
23475 d691ed308237 maint: Clean up #includes in liboctave/numeric directory. Rik <rik@octave.org> parents: 23450 diff changeset	36 #include "fDiagMatrix.h"
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	37 #include "fMatrix.h"
cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	38 #include "lo-error.h"
22322 93b3cdd36854 move most f77 function decls to separate header files John W. Eaton <jwe@octave.org> parents: 22317 diff changeset	39 #include "lo-lapack-proto.h"
23475 d691ed308237 maint: Clean up #includes in liboctave/numeric directory. Rik <rik@octave.org> parents: 23450 diff changeset	40 #include "svd.h"
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	41
22317 6ca3acf5fad8 move some new numeric classes to namespace octave::math John W. Eaton <jwe@octave.org> parents: 22297 diff changeset	42 namespace octave
6ca3acf5fad8 move some new numeric classes to namespace octave::math John W. Eaton <jwe@octave.org> parents: 22297 diff changeset	43 {
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	44 namespace math
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	45 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	46 template <typename T>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	47 T
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	48 svd<T>::left_singular_matrix (void) const
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	49 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	50 if (m_type == svd::Type::sigma_only)
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	51 (*current_liboctave_error_handler)
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	52 ("svd: U not computed because type == svd::sigma_only");
21136 7cac4e7458f2 maint: clean up code around calls to current_liboctave_error_handler. Rik <rik@octave.org> parents: 19697 diff changeset	53
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	54 return left_sm;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	55 }
1543 d6e96e0bc681 [project @ 1995-10-06 05:49:21 by jwe] jwe parents: 1368 diff changeset	56
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	57 template <typename T>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	58 T
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	59 svd<T>::right_singular_matrix (void) const
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	60 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	61 if (m_type == svd::Type::sigma_only)
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	62 (*current_liboctave_error_handler)
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	63 ("svd: V not computed because type == svd::sigma_only");
21136 7cac4e7458f2 maint: clean up code around calls to current_liboctave_error_handler. Rik <rik@octave.org> parents: 19697 diff changeset	64
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	65 return right_sm;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	66 }
1543 d6e96e0bc681 [project @ 1995-10-06 05:49:21 by jwe] jwe parents: 1368 diff changeset	67
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	68 // GESVD specializations
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	69
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	70 #define GESVD_REAL_STEP(f, F) \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	71 F77_XFCN (f, F, (F77_CONST_CHAR_ARG2 (&jobu, 1), \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	72 F77_CONST_CHAR_ARG2 (&jobv, 1), \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	73 m, n, tmp_data, m1, s_vec, u, m1, vt, \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	74 nrow_vt1, work.data (), lwork, info \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	75 F77_CHAR_ARG_LEN (1) \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	76 F77_CHAR_ARG_LEN (1)))
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	77
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	78 #define GESVD_COMPLEX_STEP(f, F, CMPLX_ARG) \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	79 F77_XFCN (f, F, (F77_CONST_CHAR_ARG2 (&jobu, 1), \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	80 F77_CONST_CHAR_ARG2 (&jobv, 1), \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	81 m, n, CMPLX_ARG (tmp_data), \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	82 m1, s_vec, CMPLX_ARG (u), m1, \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	83 CMPLX_ARG (vt), nrow_vt1, \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	84 CMPLX_ARG (work.data ()), \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	85 lwork, rwork.data (), info \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	86 F77_CHAR_ARG_LEN (1) \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	87 F77_CHAR_ARG_LEN (1)))
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	88
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	89 // DGESVD
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	90 template<>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	91 void
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	92 svd<Matrix>::gesvd (char& jobu, char& jobv, F77_INT m, F77_INT n,
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	93 double tmp_data, F77_INT m1, double s_vec,
c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	94 double u, double vt, F77_INT nrow_vt1,
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	95 std::vector<double>& work, F77_INT& lwork,
5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	96 F77_INT& info)
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	97 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	98 GESVD_REAL_STEP (dgesvd, DGESVD);
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	99
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	100 lwork = static_cast<F77_INT> (work[0]);
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	101 work.reserve (lwork);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	102
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	103 GESVD_REAL_STEP (dgesvd, DGESVD);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	104 }
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	105
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	106 // SGESVD
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	107 template<>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	108 void
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	109 svd<FloatMatrix>::gesvd (char& jobu, char& jobv, F77_INT m, F77_INT n,
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	110 float tmp_data, F77_INT m1, float s_vec,
c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	111 float u, float vt, F77_INT nrow_vt1,
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	112 std::vector<float>& work, F77_INT& lwork,
5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	113 F77_INT& info)
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	114 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	115 GESVD_REAL_STEP (sgesvd, SGESVD);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	116
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	117 lwork = static_cast<F77_INT> (work[0]);
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	118 work.reserve (lwork);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	119
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	120 GESVD_REAL_STEP (sgesvd, SGESVD);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	121 }
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	122
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	123 // ZGESVD
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	124 template<>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	125 void
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	126 svd<ComplexMatrix>::gesvd (char& jobu, char& jobv, F77_INT m, F77_INT n,
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	127 Complex tmp_data, F77_INT m1, double s_vec,
c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	128 Complex u, Complex vt, F77_INT nrow_vt1,
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	129 std::vector<Complex>& work, F77_INT& lwork,
5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	130 F77_INT& info)
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	131 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	132 std::vector<double> rwork (5 * std::max (m, n));
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	133
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	134 GESVD_COMPLEX_STEP (zgesvd, ZGESVD, F77_DBLE_CMPLX_ARG);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	135
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	136 lwork = static_cast<F77_INT> (work[0].real ());
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	137 work.reserve (lwork);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	138
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	139 GESVD_COMPLEX_STEP (zgesvd, ZGESVD, F77_DBLE_CMPLX_ARG);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	140 }
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	141
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	142 // CGESVD
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	143 template<>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	144 void
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	145 svd<FloatComplexMatrix>::gesvd (char& jobu, char& jobv, F77_INT m,
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	146 F77_INT n, FloatComplex *tmp_data,
c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	147 F77_INT m1, float s_vec, FloatComplex u,
c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	148 FloatComplex *vt, F77_INT nrow_vt1,
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	149 std::vector<FloatComplex>& work,
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	150 F77_INT& lwork, F77_INT& info)
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	151 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	152 std::vector<float> rwork (5 * std::max (m, n));
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	153
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	154 GESVD_COMPLEX_STEP (cgesvd, CGESVD, F77_CMPLX_ARG);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	155
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	156 lwork = static_cast<F77_INT> (work[0].real ());
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	157 work.reserve (lwork);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	158
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	159 GESVD_COMPLEX_STEP (cgesvd, CGESVD, F77_CMPLX_ARG);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	160 }
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	161
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	162 #undef GESVD_REAL_STEP
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	163 #undef GESVD_COMPLEX_STEP
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	164
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	165 // GESDD specializations
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	166
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	167 #define GESDD_REAL_STEP(f, F) \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	168 F77_XFCN (f, F, (F77_CONST_CHAR_ARG2 (&jobz, 1), \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	169 m, n, tmp_data, m1, s_vec, u, m1, vt, nrow_vt1, \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	170 work.data (), lwork, iwork, info \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	171 F77_CHAR_ARG_LEN (1)))
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	172
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	173 #define GESDD_COMPLEX_STEP(f, F, CMPLX_ARG) \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	174 F77_XFCN (f, F, (F77_CONST_CHAR_ARG2 (&jobz, 1), m, n, \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	175 CMPLX_ARG (tmp_data), m1, \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	176 s_vec, CMPLX_ARG (u), m1, \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	177 CMPLX_ARG (vt), nrow_vt1, \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	178 CMPLX_ARG (work.data ()), lwork, \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	179 rwork.data (), iwork, info \
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	180 F77_CHAR_ARG_LEN (1)))
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	181
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	182 // DGESDD
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	183 template<>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	184 void
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	185 svd<Matrix>::gesdd (char& jobz, F77_INT m, F77_INT n, double *tmp_data,
c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	186 F77_INT m1, double s_vec, double u, double *vt,
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	187 F77_INT nrow_vt1, std::vector<double>& work,
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	188 F77_INT& lwork, F77_INT *iwork, F77_INT& info)
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	189 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	190 GESDD_REAL_STEP (dgesdd, DGESDD);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	191
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	192 lwork = static_cast<F77_INT> (work[0]);
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	193 work.reserve (lwork);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	194
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	195 GESDD_REAL_STEP (dgesdd, DGESDD);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	196 }
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	197
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	198 // SGESDD
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	199 template<>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	200 void
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	201 svd<FloatMatrix>::gesdd (char& jobz, F77_INT m, F77_INT n, float *tmp_data,
c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	202 F77_INT m1, float s_vec, float u, float *vt,
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	203 F77_INT nrow_vt1, std::vector<float>& work,
5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	204 F77_INT& lwork, F77_INT* iwork, F77_INT& info)
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	205 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	206 GESDD_REAL_STEP (sgesdd, SGESDD);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	207
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	208 lwork = static_cast<F77_INT> (work[0]);
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	209 work.reserve (lwork);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	210
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	211 GESDD_REAL_STEP (sgesdd, SGESDD);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	212 }
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	213
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	214 // ZGESDD
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	215 template<>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	216 void
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	217 svd<ComplexMatrix>::gesdd (char& jobz, F77_INT m, F77_INT n,
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	218 Complex tmp_data, F77_INT m1, double s_vec,
c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	219 Complex u, Complex vt, F77_INT nrow_vt1,
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	220 std::vector<Complex>& work, F77_INT& lwork,
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	221 F77_INT *iwork, F77_INT& info)
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	222 {
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	223
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	224 F77_INT min_mn = std::min (m, n);
25349 a8a4b6e6e754 Fix lrwork dimension in svd for complex cases (bug #53823). Marco Caliari <marco.caliari@univr.it> parents: 25054 diff changeset	225 F77_INT max_mn = std::max (m, n);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	226
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	227 F77_INT lrwork;
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	228 if (jobz == 'N')
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	229 lrwork = 7*min_mn;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	230 else
25349 a8a4b6e6e754 Fix lrwork dimension in svd for complex cases (bug #53823). Marco Caliari <marco.caliari@univr.it> parents: 25054 diff changeset	231 lrwork = min_mn * std::max (5min_mn+5, 2max_mn+2*min_mn+1);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	232
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	233 std::vector<double> rwork (lrwork);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	234
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	235 GESDD_COMPLEX_STEP (zgesdd, ZGESDD, F77_DBLE_CMPLX_ARG);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	236
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	237 lwork = static_cast<F77_INT> (work[0].real ());
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	238 work.reserve (lwork);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	239
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	240 GESDD_COMPLEX_STEP (zgesdd, ZGESDD, F77_DBLE_CMPLX_ARG);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	241 }
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	242
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	243 // CGESDD
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	244 template<>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	245 void
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	246 svd<FloatComplexMatrix>::gesdd (char& jobz, F77_INT m, F77_INT n,
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	247 FloatComplex *tmp_data, F77_INT m1,
c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	248 float s_vec, FloatComplex u,
c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	249 FloatComplex *vt, F77_INT nrow_vt1,
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	250 std::vector<FloatComplex>& work,
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	251 F77_INT& lwork, F77_INT *iwork,
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	252 F77_INT& info)
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	253 {
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	254 F77_INT min_mn = std::min (m, n);
5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	255 F77_INT max_mn = std::max (m, n);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	256
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	257 F77_INT lrwork;
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	258 if (jobz == 'N')
25349 a8a4b6e6e754 Fix lrwork dimension in svd for complex cases (bug #53823). Marco Caliari <marco.caliari@univr.it> parents: 25054 diff changeset	259 lrwork = 7*min_mn;
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	260 else
25349 a8a4b6e6e754 Fix lrwork dimension in svd for complex cases (bug #53823). Marco Caliari <marco.caliari@univr.it> parents: 25054 diff changeset	261 lrwork = min_mn * std::max (5min_mn+5, 2max_mn+2*min_mn+1);
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	262 std::vector<float> rwork (lrwork);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	263
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	264 GESDD_COMPLEX_STEP (cgesdd, CGESDD, F77_CMPLX_ARG);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	265
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	266 lwork = static_cast<F77_INT> (work[0].real ());
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	267 work.reserve (lwork);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	268
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	269 GESDD_COMPLEX_STEP (cgesdd, CGESDD, F77_CMPLX_ARG);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	270 }
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	271
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	272 #undef GESDD_REAL_STEP
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	273 #undef GESDD_COMPLEX_STEP
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	274
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	275 template<typename T>
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	276 svd<T>::svd (const T& a, svd::Type type,
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	277 svd::Driver driver)
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	278 : m_type (type), m_driver (driver), left_sm (), sigma (), right_sm ()
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	279 {
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	280 F77_INT info;
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	281
23696 08036a7f3660 remove octave:: namespace tag from symbols used inside octave namespace John W. Eaton <jwe@octave.org> parents: 23475 diff changeset	282 F77_INT m = to_f77_int (a.rows ());
08036a7f3660 remove octave:: namespace tag from symbols used inside octave namespace John W. Eaton <jwe@octave.org> parents: 23475 diff changeset	283 F77_INT n = to_f77_int (a.cols ());
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	284
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	285 if (m == 0 \|\| n == 0)
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	286 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	287 switch (m_type)
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	288 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	289 case svd::Type::std:
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	290 left_sm = T (m, m, 0);
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	291 for (F77_INT i = 0; i < m; i++)
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	292 left_sm.xelem (i, i) = 1;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	293 sigma = DM_T (m, n);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	294 right_sm = T (n, n, 0);
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	295 for (F77_INT i = 0; i < n; i++)
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	296 right_sm.xelem (i, i) = 1;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	297 break;
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	298
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	299 case svd::Type::economy:
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	300 left_sm = T (m, 0, 0);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	301 sigma = DM_T (0, 0);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	302 right_sm = T (0, n, 0);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	303 break;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	304
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	305 case svd::Type::sigma_only:
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	306 default:
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	307 sigma = DM_T (0, 1);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	308 break;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	309 }
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	310 return;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	311 }
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	312
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	313 T atmp = a;
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	314 P *tmp_data = atmp.fortran_vec ();
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	315
23450 855122b993da maint: Wrap tertiary operator in parentheses "(COND ? x : y)". Rik <rik@octave.org> parents: 23449 diff changeset	316 F77_INT min_mn = (m < n ? m : n);
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	317
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	318 char jobu = 'A';
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	319 char jobv = 'A';
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	320
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	321 F77_INT ncol_u = m;
5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	322 F77_INT nrow_vt = n;
5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	323 F77_INT nrow_s = m;
5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	324 F77_INT ncol_s = n;
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	325
22228 4afe3705ea75 svd: do not shadow member variables with local variables. Carnë Draug <carandraug@octave.org> parents: 22204 diff changeset	326 switch (m_type)
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	327 {
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	328 case svd::Type::economy:
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	329 jobu = jobv = 'S';
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	330 ncol_u = nrow_vt = nrow_s = ncol_s = min_mn;
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	331 break;
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	332
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	333 case svd::Type::sigma_only:
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	334
23475 d691ed308237 maint: Clean up #includes in liboctave/numeric directory. Rik <rik@octave.org> parents: 23450 diff changeset	335 // Note: for this case, both jobu and jobv should be 'N', but there
d691ed308237 maint: Clean up #includes in liboctave/numeric directory. Rik <rik@octave.org> parents: 23450 diff changeset	336 // seems to be a bug in dgesvd from Lapack V2.0. To demonstrate the
d691ed308237 maint: Clean up #includes in liboctave/numeric directory. Rik <rik@octave.org> parents: 23450 diff changeset	337 // bug, set both jobu and jobv to 'N' and find the singular values of
d691ed308237 maint: Clean up #includes in liboctave/numeric directory. Rik <rik@octave.org> parents: 23450 diff changeset	338 // [eye(3), eye(3)]. The result is [-sqrt(2), -sqrt(2), -sqrt(2)].
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	339 //
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	340 // For Lapack 3.0, this problem seems to be fixed.
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	341
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	342 jobu = jobv = 'N';
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	343 ncol_u = nrow_vt = 1;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	344 break;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	345
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	346 default:
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	347 break;
469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	348 }
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	349
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	350 if (! (jobu == 'N' \|\| jobu == 'O'))
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	351 left_sm.resize (m, ncol_u);
22204 469c817eb256 svd: reduce code duplication with more use of template and macro. Carnë Draug <carandraug@octave.org> parents: 22135 diff changeset	352
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	353 P *u = left_sm.fortran_vec ();
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	354
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	355 sigma.resize (nrow_s, ncol_s);
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	356 DM_P *s_vec = sigma.fortran_vec ();
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	357
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	358 if (! (jobv == 'N' \|\| jobv == 'O'))
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	359 right_sm.resize (nrow_vt, n);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	360
23449 c763214a8260 maint: Use convention 'int x' for naming pointers. Rik <rik@octave.org>* parents: 23448 diff changeset	361 P *vt = right_sm.fortran_vec ();
537 4ecbfd3c3710 [project @ 1994-07-21 22:30:34 by jwe] jwe parents: 457 diff changeset	362
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	363 // Query _GESVD for the correct dimension of WORK.
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	364
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	365 F77_INT lwork = -1;
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	366
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	367 std::vector<P> work (1);
537 4ecbfd3c3710 [project @ 1994-07-21 22:30:34 by jwe] jwe parents: 457 diff changeset	368
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	369 F77_INT m1 = std::max (m, static_cast<F77_INT> (1));
5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	370 F77_INT nrow_vt1 = std::max (nrow_vt, static_cast<F77_INT> (1));
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	371
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	372 if (m_driver == svd::Driver::GESVD)
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	373 gesvd (jobu, jobv, m, n, tmp_data, m1, s_vec, u, vt, nrow_vt1,
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	374 work, lwork, info);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	375 else if (m_driver == svd::Driver::GESDD)
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	376 {
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	377 assert (jobu == jobv);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	378 char jobz = jobu;
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	379
22944 5126040e8f49 use F77_INT instead of octave_idx_type for liboctave svd class John W. Eaton <jwe@octave.org> parents: 22755 diff changeset	380 std::vector<F77_INT> iwork (8 * std::min (m, n));
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	381
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	382 gesdd (jobz, m, n, tmp_data, m1, s_vec, u, vt, nrow_vt1,
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	383 work, lwork, iwork.data (), info);
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	384 }
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	385 else
23448 e1c02e7126a9 avoid calling abort if possible; include cstdlib for declaration if needed John W. Eaton <jwe@octave.org> parents: 23220 diff changeset	386 (*current_liboctave_error_handler) ("svd: unknown driver");
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	387
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	388 if (! (jobv == 'N' \|\| jobv == 'O'))
22592 ee0df00e12d6 fix complex svd (bug #49232) John W. Eaton <jwe@octave.org> parents: 22402 diff changeset	389 right_sm = right_sm.hermitian ();
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	390 }
cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	391
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	392 // Instantiations we need.
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	393
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	394 template class svd<Matrix>;
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	395
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	396 template class svd<FloatMatrix>;
21273 cbced1c09916 better use of templates for svd classes John W. Eaton <jwe@octave.org> parents: 21202 diff changeset	397
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	398 template class svd<ComplexMatrix>;
3336 08ad797989f8 [project @ 1999-11-03 21:41:34 by jwe] jwe parents: 3335 diff changeset	399
22329 7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	400 template class svd<FloatComplexMatrix>;
7f3c7a8bd131 maint: Indent namespaces in liboctave/numeric files. John W. Eaton <jwe@octave.org> parents: 22323 diff changeset	401 }
457 3d4b4f0fa5ba [project @ 1994-06-06 00:33:33 by jwe] jwe parents: diff changeset	402 }

Mercurial > octave

annotate liboctave/numeric/svd.cc @ 25349:a8a4b6e6e754 stable