home *** CD-ROM | disk | FTP | other *** search
/ APDL Eductation Resources / APDL Eductation Resources.iso / programs / electronic / rlab / TestMatrix / kron_r < prev    next >
Encoding:
Text File  |  1994-12-20  |  1.6 KB  |  67 lines
  1. //-------------------------------------------------------------------//
  2.  
  3. // Synopsis:    Kronecker tensor product.
  4.  
  5. // Syntax:    kron ( x , y )
  6.  
  7. // Description:
  8.  
  9. //    kron(X,Y) is the Kronecker tensor product of X and Y. The
  10. //    result is a large matrix formed by taking all possible
  11. //    products between the elements of X and those of Y.   For
  12. //    example, if X is 2 by 3, then KRON(X,Y) is 
  13. //
  14. //       [ X[1;1]*Y  X[1;2]*Y  X[1;3]*Y
  15. //         X[2;1]*Y  X[2;2]*Y  X[2;3]*Y ]
  16. //
  17. //    If either X or Y is sparse, only nonzero elements are
  18. //    multiplied in the computation, and the result is sparse. 
  19.  
  20. //    Full: J. N. Little, 4-21-85.
  21. //    Sparse: T. R. Gardos (Georgia Tech), 5-4-93.
  22. //    Copyright (c) 1984-93 by The MathWorks, Inc.
  23. //
  24. //      This version not strictly part of the toolbox - included here
  25. //      because version supplied with Matlab 4.0 is not sparse aware.
  26. //      (This comment does not apply to RLaB)
  27.  
  28. //    This file is a translation of kron.m from version 2.0 of
  29. //    "The Test Matrix Toolbox for Matlab", described in Numerical
  30. //    Analysis Report No. 237, December 1993, by N. J. Higham.
  31.  
  32. //-------------------------------------------------------------------//
  33.  
  34. kron = function ( A , B )
  35. {
  36.   ma = A.nr; na = A.nr;
  37.   mb = A.nr; nb = B.nc;
  38.  
  39.   // full inputs
  40.  
  41.   K = zeros(ma*mb,na*nb);
  42.   if (ma*na <= mb*nb)
  43.   {
  44.     for (i in 1:ma)
  45.     {
  46.       ik = 1+(i-1)*mb:i*mb;
  47.       for (j in 1:na)
  48.       {
  49.         jk = 1+(j-1)*nb:j*nb;
  50.         K[ik;jk] = A[i;j]*B;
  51.       }
  52.     }
  53.   else
  54.     for (i in 1:mb)
  55.     {
  56.       ik = i:mb:(ma-1)*mb+i;
  57.       for (j in 1:nb)
  58.       {
  59.         jk = j:nb:(na-1)*nb+j;
  60.         K[ik;jk] = B[i;j]*A;
  61.       }
  62.     }
  63.   }
  64.  
  65.   return K;
  66. };
  67.