# HG changeset patch
# User Mike Miller <mtmiller@octave.org>
# Date 1442321315 14400
# Node ID ff904ae0285b5415947e4b127fd8bffab7b8dc4a
# Parent  ba4088aee342040e7a9148265ffff23cba4259a6
eig: Return correct solution for a pair of hermitian matrices (bug #45511)

* liboctave/numeric/EIG.cc (EIG::init): Use correct form of hermitian_init.
* liboctave/numeric/fEIG.cc (FloatEIG::init): Likewise.
* libinterp/corefcn/eig.cc: Add %!test cases.

Thanks to Marco Caliari for identifying the fix.

diff -r ba4088aee342 -r ff904ae0285b libinterp/corefcn/eig.cc
--- a/libinterp/corefcn/eig.cc	Sun Sep 13 14:47:58 2015 -0700
+++ b/libinterp/corefcn/eig.cc	Tue Sep 15 08:48:35 2015 -0400
@@ -324,6 +324,18 @@
 %! assert (A * v(:, 1), d(1, 1) * B * v(:, 1), sqrt (eps));
 %! assert (A * v(:, 2), d(2, 2) * B * v(:, 2), sqrt (eps));
 
+%!test
+%! A = [1, 1+i; 1-i, 1];  B = [2, 0; 0, 2];
+%! [v, d] = eig (A, B);
+%! assert (A * v(:, 1), d(1, 1) * B * v(:, 1), sqrt (eps));
+%! assert (A * v(:, 2), d(2, 2) * B * v(:, 2), sqrt (eps));
+
+%!test
+%! A = single ([1, 1+i; 1-i, 1]);  B = single ([2, 0; 0, 2]);
+%! [v, d] = eig (A, B);
+%! assert (A * v(:, 1), d(1, 1) * B * v(:, 1), sqrt (eps ("single")));
+%! assert (A * v(:, 2), d(2, 2) * B * v(:, 2), sqrt (eps ("single")));
+
 %!error eig ()
 %!error eig ([1, 2; 3, 4], [4, 3; 2, 1], 1)
 %!error <EIG requires same size matrices> eig ([1, 2; 3, 4], 2)
diff -r ba4088aee342 -r ff904ae0285b liboctave/numeric/EIG.cc
--- a/liboctave/numeric/EIG.cc	Sun Sep 13 14:47:58 2015 -0700
+++ b/liboctave/numeric/EIG.cc	Tue Sep 15 08:48:35 2015 -0400
@@ -715,7 +715,7 @@
                              F77_CHAR_ARG_LEN (1)));
 
   if (a.is_hermitian () && b.is_hermitian () && info == 0)
-    return hermitian_init (a, calc_ev);
+    return hermitian_init (a, b, calc_ev);
 
   ComplexMatrix atmp = a;
   Complex *atmp_data = atmp.fortran_vec ();
diff -r ba4088aee342 -r ff904ae0285b liboctave/numeric/fEIG.cc
--- a/liboctave/numeric/fEIG.cc	Sun Sep 13 14:47:58 2015 -0700
+++ b/liboctave/numeric/fEIG.cc	Tue Sep 15 08:48:35 2015 -0400
@@ -712,7 +712,7 @@
                              F77_CHAR_ARG_LEN (1)));
 
   if (a.is_hermitian () && b.is_hermitian () && info == 0)
-    return hermitian_init (a, calc_ev);
+    return hermitian_init (a, b, calc_ev);
 
   FloatComplexMatrix atmp = a;
   FloatComplex *atmp_data = atmp.fortran_vec ();