# HG changeset patch # User John W. Eaton # Date 1205246973 14400 # Node ID b7a30502f0c965e8988dd85d4a7da29d4260d7cd # Parent 114dcfc154325a737c77581c1b38ba1401667f06 handle possible error from EIG diff --git a/src/ChangeLog b/src/ChangeLog --- a/src/ChangeLog +++ b/src/ChangeLog @@ -1,3 +1,9 @@ +2008-03-11 John W. Eaton + + * DLD-FUNCTIONS/eig.cc (Feig): Handle possible error from EIG. + * DLD-FUNCTIONS/qp.cc (qp, Fqp): Likewise. + * xpow.cc (xpow): Likewise. + 2008-03-10 John W. Eaton * mex.cc (mxCreateLogicalScalar): Argument is now mxLogical. diff --git a/src/DLD-FUNCTIONS/__qp__.cc b/src/DLD-FUNCTIONS/__qp__.cc --- a/src/DLD-FUNCTIONS/__qp__.cc +++ b/src/DLD-FUNCTIONS/__qp__.cc @@ -137,6 +137,13 @@ // Computing the ??? EIG eigH (H); + + if (error_state) + { + error ("qp: failed to compute eigenvalues of H"); + return -1; + } + ColumnVector eigenvalH = real (eigH.eigenvalues ()); Matrix eigenvecH = real (eigH.eigenvectors ()); double minReal = eigenvalH.min (); @@ -272,6 +279,13 @@ // Searching for the most negative curvature. EIG eigrH (rH); + + if (error_state) + { + error ("qp: failed to compute eigenvalues of rH"); + return -1; + } + ColumnVector eigenvalrH = real (eigrH.eigenvalues ()); Matrix eigenvecrH = real (eigrH.eigenvectors ()); double mRrH = eigenvalrH.min (); @@ -494,10 +508,15 @@ int info = qp (H, q, Aeq, beq, Ain, bin, maxit, x, lambda, iter); - retval(3) = iter; - retval(2) = info; - retval(1) = lambda; - retval(0) = x; + if (! error_state) + { + retval(3) = iter; + retval(2) = info; + retval(1) = lambda; + retval(0) = x; + } + else + error ("qp: internal error"); } else error ("__qp__: invalid arguments"); diff --git a/src/DLD-FUNCTIONS/eig.cc b/src/DLD-FUNCTIONS/eig.cc --- a/src/DLD-FUNCTIONS/eig.cc +++ b/src/DLD-FUNCTIONS/eig.cc @@ -101,18 +101,21 @@ return retval; } - if (nargout == 0 || nargout == 1) - { - retval(0) = result.eigenvalues (); - } - else + if (! error_state) { - // Blame it on Matlab. + if (nargout == 0 || nargout == 1) + { + retval(0) = result.eigenvalues (); + } + else + { + // Blame it on Matlab. - ComplexDiagMatrix d (result.eigenvalues ()); + ComplexDiagMatrix d (result.eigenvalues ()); - retval(1) = d; - retval(0) = result.eigenvectors (); + retval(1) = d; + retval(0) = result.eigenvectors (); + } } return retval; diff --git a/src/xpow.cc b/src/xpow.cc --- a/src/xpow.cc +++ b/src/xpow.cc @@ -101,20 +101,26 @@ else { EIG b_eig (b); - ComplexColumnVector lambda (b_eig.eigenvalues ()); - ComplexMatrix Q (b_eig.eigenvectors ()); - for (octave_idx_type i = 0; i < nr; i++) + if (! error_state) { - Complex elt = lambda (i); - if (std::imag (elt) == 0.0) - lambda (i) = std::pow (a, std::real (elt)); - else - lambda (i) = std::pow (a, elt); + ComplexColumnVector lambda (b_eig.eigenvalues ()); + ComplexMatrix Q (b_eig.eigenvectors ()); + + for (octave_idx_type i = 0; i < nr; i++) + { + Complex elt = lambda(i); + if (std::imag (elt) == 0.0) + lambda(i) = std::pow (a, std::real (elt)); + else + lambda(i) = std::pow (a, elt); + } + ComplexDiagMatrix D (lambda); + + retval = ComplexMatrix (Q * D * Q.inverse ()); } - ComplexDiagMatrix D (lambda); - - retval = ComplexMatrix (Q * D * Q.inverse ()); + else + error ("xpow: matrix diagonalization failed"); } return retval; @@ -144,20 +150,26 @@ else { EIG b_eig (b); - ComplexColumnVector lambda (b_eig.eigenvalues ()); - ComplexMatrix Q (b_eig.eigenvectors ()); - for (octave_idx_type i = 0; i < nr; i++) + if (! error_state) { - Complex elt = lambda (i); - if (std::imag (elt) == 0.0) - lambda (i) = std::pow (a, std::real (elt)); - else - lambda (i) = std::pow (a, elt); + ComplexColumnVector lambda (b_eig.eigenvalues ()); + ComplexMatrix Q (b_eig.eigenvectors ()); + + for (octave_idx_type i = 0; i < nr; i++) + { + Complex elt = lambda(i); + if (std::imag (elt) == 0.0) + lambda(i) = std::pow (a, std::real (elt)); + else + lambda(i) = std::pow (a, elt); + } + ComplexDiagMatrix D (lambda); + + retval = ComplexMatrix (Q * D * Q.inverse ()); } - ComplexDiagMatrix D (lambda); - - retval = ComplexMatrix (Q * D * Q.inverse ()); + else + error ("xpow: matrix diagonalization failed"); } return retval; @@ -228,15 +240,21 @@ else { EIG a_eig (a); - ComplexColumnVector lambda (a_eig.eigenvalues ()); - ComplexMatrix Q (a_eig.eigenvectors ()); + + if (! error_state) + { + ComplexColumnVector lambda (a_eig.eigenvalues ()); + ComplexMatrix Q (a_eig.eigenvectors ()); - for (octave_idx_type i = 0; i < nr; i++) - lambda (i) = std::pow (lambda (i), b); + for (octave_idx_type i = 0; i < nr; i++) + lambda(i) = std::pow (lambda(i), b); + + ComplexDiagMatrix D (lambda); - ComplexDiagMatrix D (lambda); - - retval = ComplexMatrix (Q * D * Q.inverse ()); + retval = ComplexMatrix (Q * D * Q.inverse ()); + } + else + error ("xpow: matrix diagonalization failed"); } } @@ -257,15 +275,21 @@ else { EIG a_eig (a); - ComplexColumnVector lambda (a_eig.eigenvalues ()); - ComplexMatrix Q (a_eig.eigenvectors ()); + + if (! error_state) + { + ComplexColumnVector lambda (a_eig.eigenvalues ()); + ComplexMatrix Q (a_eig.eigenvectors ()); - for (octave_idx_type i = 0; i < nr; i++) - lambda (i) = std::pow (lambda (i), b); + for (octave_idx_type i = 0; i < nr; i++) + lambda(i) = std::pow (lambda(i), b); + + ComplexDiagMatrix D (lambda); - ComplexDiagMatrix D (lambda); - - retval = ComplexMatrix (Q * D * Q.inverse ()); + retval = ComplexMatrix (Q * D * Q.inverse ()); + } + else + error ("xpow: matrix diagonalization failed"); } return retval; @@ -299,20 +323,26 @@ else { EIG b_eig (b); - ComplexColumnVector lambda (b_eig.eigenvalues ()); - ComplexMatrix Q (b_eig.eigenvectors ()); - for (octave_idx_type i = 0; i < nr; i++) + if (! error_state) { - Complex elt = lambda (i); - if (std::imag (elt) == 0.0) - lambda (i) = std::pow (a, std::real (elt)); - else - lambda (i) = std::pow (a, elt); + ComplexColumnVector lambda (b_eig.eigenvalues ()); + ComplexMatrix Q (b_eig.eigenvectors ()); + + for (octave_idx_type i = 0; i < nr; i++) + { + Complex elt = lambda(i); + if (std::imag (elt) == 0.0) + lambda(i) = std::pow (a, std::real (elt)); + else + lambda(i) = std::pow (a, elt); + } + ComplexDiagMatrix D (lambda); + + retval = ComplexMatrix (Q * D * Q.inverse ()); } - ComplexDiagMatrix D (lambda); - - retval = ComplexMatrix (Q * D * Q.inverse ()); + else + error ("xpow: matrix diagonalization failed"); } return retval; @@ -341,20 +371,26 @@ else { EIG b_eig (b); - ComplexColumnVector lambda (b_eig.eigenvalues ()); - ComplexMatrix Q (b_eig.eigenvectors ()); - for (octave_idx_type i = 0; i < nr; i++) + if (! error_state) { - Complex elt = lambda (i); - if (std::imag (elt) == 0.0) - lambda (i) = std::pow (a, std::real (elt)); - else - lambda (i) = std::pow (a, elt); + ComplexColumnVector lambda (b_eig.eigenvalues ()); + ComplexMatrix Q (b_eig.eigenvectors ()); + + for (octave_idx_type i = 0; i < nr; i++) + { + Complex elt = lambda(i); + if (std::imag (elt) == 0.0) + lambda(i) = std::pow (a, std::real (elt)); + else + lambda(i) = std::pow (a, elt); + } + ComplexDiagMatrix D (lambda); + + retval = ComplexMatrix (Q * D * Q.inverse ()); } - ComplexDiagMatrix D (lambda); - - retval = ComplexMatrix (Q * D * Q.inverse ()); + else + error ("xpow: matrix diagonalization failed"); } return retval; @@ -425,15 +461,21 @@ else { EIG a_eig (a); - ComplexColumnVector lambda (a_eig.eigenvalues ()); - ComplexMatrix Q (a_eig.eigenvectors ()); + + if (! error_state) + { + ComplexColumnVector lambda (a_eig.eigenvalues ()); + ComplexMatrix Q (a_eig.eigenvectors ()); - for (octave_idx_type i = 0; i < nr; i++) - lambda (i) = std::pow (lambda (i), b); + for (octave_idx_type i = 0; i < nr; i++) + lambda(i) = std::pow (lambda(i), b); + + ComplexDiagMatrix D (lambda); - ComplexDiagMatrix D (lambda); - - retval = ComplexMatrix (Q * D * Q.inverse ()); + retval = ComplexMatrix (Q * D * Q.inverse ()); + } + else + error ("xpow: matrix diagonalization failed"); } } @@ -454,15 +496,21 @@ else { EIG a_eig (a); - ComplexColumnVector lambda (a_eig.eigenvalues ()); - ComplexMatrix Q (a_eig.eigenvectors ()); + + if (! error_state) + { + ComplexColumnVector lambda (a_eig.eigenvalues ()); + ComplexMatrix Q (a_eig.eigenvectors ()); - for (octave_idx_type i = 0; i < nr; i++) - lambda (i) = std::pow (lambda (i), b); + for (octave_idx_type i = 0; i < nr; i++) + lambda(i) = std::pow (lambda(i), b); + + ComplexDiagMatrix D (lambda); - ComplexDiagMatrix D (lambda); - - retval = ComplexMatrix (Q * D * Q.inverse ()); + retval = ComplexMatrix (Q * D * Q.inverse ()); + } + else + error ("xpow: matrix diagonalization failed"); } return retval;