Mercurial > hg > octave-lyh

changeset 11682:b7a30502f0c9 release-3-0-x

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
handle possible error from EIG
author John W. Eaton <jwe@octave.org>
date Tue, 11 Mar 2008 10:49:33 -0400
parents 114dcfc15432
children 738ca2f7c83e
files src/ChangeLog src/DLD-FUNCTIONS/__qp__.cc src/DLD-FUNCTIONS/eig.cc src/xpow.cc
diffstat 4 files changed, 161 insertions(+), 85 deletions(-) [+]
line wrap: on
line diff
--- a/src/ChangeLog
+++ b/src/ChangeLog
@@ -1,3 +1,9 @@
+2008-03-11  John W. Eaton  <jwe@octave.org>
+
+	* 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  <jwe@octave.org>
 
 	* mex.cc (mxCreateLogicalScalar): Argument is now mxLogical.
--- 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");
--- 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;
--- 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;