changeset 14555:6eabd81604b5 stable

allow kron to work for two diag matrix arguments (bug #35647) * kron.cc (dispatch_kron): Fix recursive call for case of two diagonal matrix objects as arguments. New tests.
author John W. Eaton <jwe@octave.org>
date Thu, 12 Apr 2012 14:01:39 -0400
parents 93cb513ed5fc
children 15e4ec503cfd 05b7fa3b64c8
files src/DLD-FUNCTIONS/kron.cc
diffstat 1 files changed, 25 insertions(+), 8 deletions(-) [+]
line wrap: on
line diff
--- a/src/DLD-FUNCTIONS/kron.cc	Thu Apr 05 22:29:28 2012 -0700
+++ b/src/DLD-FUNCTIONS/kron.cc	Thu Apr 12 14:01:39 2012 -0400
@@ -162,7 +162,6 @@
   return c;
 }
 
-
 template <class MTA, class MTB>
 octave_value
 do_kron (const octave_value& a, const octave_value& b)
@@ -183,12 +182,18 @@
       if (b.is_diag_matrix () && a.rows () == a.columns ()
           && b.rows () == b.columns ())
         {
-          octave_value_list tmp;
-          tmp(0) = a.diag ();
-          tmp(1) = b.diag ();
-          tmp = dispatch_kron (tmp, 1);
-          if (tmp.length () == 1)
-            retval = tmp(0).diag ();
+          // We have two diagonal matrices, the product of those will be
+          // another diagonal matrix.  To do that efficiently, extract
+          // the diagonals as vectors and compute the product.  That
+          // will be another vector, which we then use to construct a
+          // diagonal matrix object.  Note that this will fail if our
+          // digaonal matrix object is modified to allow the non-zero
+          // values to be stored off of the principal diagonal (i.e., if
+          // diag ([1,2], 3) is modified to return a diagonal matrix
+          // object instead of a full matrix object).
+
+          octave_value tmp = dispatch_kron (a.diag (), b.diag ());
+          retval = tmp.diag ();
         }
       else if (a.is_single_type () || b.is_single_type ())
         {
@@ -289,7 +294,6 @@
 
 
 /*
-
 %!test
 %! x = ones(2);
 %! assert( kron (x, x), ones (4));
@@ -302,4 +306,17 @@
 %!assert (kron (x, y, z), kron (kron (x, y), z))
 %!assert (kron (x, y, z), kron (x, kron (y, z)))
 
+
+%!assert (kron (diag ([1, 2]), diag ([3, 4])), diag ([3, 4, 6, 8]))
+
+%% Test for two diag matrices.  See the comments above in
+%% dispatch_kron for this case.
+%%
+%!test
+%! expected = zeros (16, 16);
+%! expected (1, 11) = 3;
+%! expected (2, 12) = 4;
+%! expected (5, 15) = 6;
+%! expected (6, 16) = 8;
+%! assert (kron (diag ([1, 2], 2), diag ([3, 4], 2)), expected)
 */