home *** CD-ROM | disk | FTP | other *** search
/ io Programmo 21 / IOPROG_21.ISO / SOFT / LIBMAT.ZIP / MATKRON.CPP < prev    next >
Encoding:
C/C++ Source or Header  |  1992-10-16  |  4.5 KB  |  144 lines
  1. /**************************************************/
  2. /*    matkron.c source for Kronecker products     */
  3. /**************************************************/
  4.  
  5.  
  6. /**************************************************/
  7. /*            MatClass Source File                */
  8. /*       Copyright of C. R. Birchenhall           */
  9. /*       University of Manchester, UK.            */
  10. /*   MatClass is freeware. This file should be    */
  11. /* made freely available to users of any software */
  12. /* whose creation is wholly or partly dependent   */
  13. /*                on this file.                   */
  14. /**************************************************/
  15.  
  16. #include "matrix.hpp"
  17.  
  18. matrix& matrix::kronOf( const matrix& x, const matrix& y,
  19.                         transposition opt )
  20. {
  21.    static char *mName = "kronOf" ;
  22.    matFunc func( mName ) ; debugInfo( func ) ;
  23.    INDEX xrow = x.nRows(), xcol = x.nCols() ;
  24.    INDEX yrow = y.nRows(), ycol = y.nCols() ;
  25.    INDEX nr, nc, i, j, zi, zj, xi, xj, yi, yj ;
  26.    REAL r ;
  27.    matError error = NOERROR ;
  28.    switch ( opt ) {
  29.       case TRANS0 :
  30.          nr = xrow * yrow ;
  31.          nc = xcol * ycol ;
  32.          break ;
  33.       case TRANS1 :
  34.          nr = xcol * yrow ;
  35.          nc = xrow * ycol ;
  36.          break ;
  37.       case TRANS2 :
  38.          nr = xrow * ycol ;
  39.          nc = xcol * yrow ;
  40.          break ;
  41.       case TRANS12 :
  42.          nr = xcol * ycol ;
  43.          nc = xrow * yrow ;
  44.          break ;
  45.       default :
  46.          error = NGPAR ;
  47.          break ;
  48.    } // switch
  49.    if ( error )
  50.      errorExit( mName, error ) ;
  51.    reset( nr, nc ) ;
  52.    if ( opt == TRANS0 ) {
  53.       for ( xi = 1, i = 1 ; xi <= xrow ; xi++, i += yrow ) {
  54.          for ( xj = 1, j = 1 ; xj <= xcol ; xj++, j += ycol ) {
  55.             r = x( xi, xj ) ;
  56.             for ( yi = 1, zi = i ; yi <= yrow ; yi++, zi++ )
  57.                for ( yj = 1, zj = j ; yj <= ycol ; yj++, zj++ )
  58.                   mat( zi, zj ) = r * y( yi, yj ) ;
  59.          } // for xj
  60.       } // for xi
  61.    } else if ( opt == TRANS1 ) {
  62.       for ( xi = 1, i = 1 ; xi <= xcol ; xi++, i += yrow ) {
  63.          for ( xj = 1, j = 1 ; xj <= xrow ; xj++, j += ycol ) {
  64.             r = x( xj, xi ) ;
  65.             for ( yi = 1, zi = i ; yi <= yrow ; yi++, zi++ )
  66.                for ( yj = 1, zj = j ; yj <= ycol ; yj++, zj++ )
  67.                   mat( zi, zj ) = r * y( yi, yj ) ;
  68.          } // for xj
  69.       } // for xi
  70.    } else if ( opt == TRANS2 ) {
  71.       for ( xi = 1, i = 1 ; xi <= xrow ; xi++, i += ycol ) {
  72.          for ( xj = 1, j = 1 ; xj <= xcol ; xj++, j += yrow ) {
  73.             r = x( xi, xj ) ;
  74.             for ( yi = 1, zi = i ; yi <= ycol ; yi++, zi++ )
  75.                for ( yj = 1, zj = j ; yj <= yrow ; yj++, zj++ )
  76.                   mat( zi, zj ) = r * y( yj, yi ) ;
  77.          } // for xj
  78.       } // for xi
  79.    } else { //  opt == TRANS12
  80.       for ( xi = 1, i = 1 ; xi <= xcol ; xi++, i += ycol ) {
  81.          for ( xj = 1 , j = 1 ; xj <= xrow ; xj++, j += yrow ) {
  82.             r = x( xj, xi ) ;
  83.             for ( yi = 1, zi = i ; yi <= ycol ; yi++, zi++ )
  84.                for ( yj = 1, zj = j ; yj <= yrow ; yj++, zj++ )
  85.                   mat( zi, zj ) = r * y( yj, yi ) ;
  86.          } // for xj
  87.       } // for xi
  88.    } // else
  89.    return *this ;
  90. } // kronOf
  91.  
  92. matrix matrix::kron( const matrix& y ) M_CONST
  93. {
  94.    static char *mName = "kron" ;
  95.    matFunc func( mName ) ; debugInfo( func ) ;
  96.    matrix z ;
  97.    z.kronOf( *this, y ) ;
  98.    return z ;
  99. } // matrix::kron
  100.  
  101. matrix matrix::TKron( const matrix& y ) M_CONST
  102. {
  103.    static char *mName = "TKron" ;
  104.    matFunc func( mName ) ; debugInfo( func ) ;
  105.    matrix z ;
  106.    z.kronOf( *this, y, TRANS1 ) ;
  107.    return z ;
  108. } // matrix::TKron
  109.  
  110. matrix matrix::kronT( const matrix& y ) M_CONST
  111. {
  112.    static char *mName = "kronT" ;
  113.    matFunc func( mName ) ; debugInfo( func ) ;
  114.    matrix z ;
  115.    z.kronOf( *this, y, TRANS2 ) ;
  116.    return z ;
  117. } // matrix::kronT
  118.  
  119. matrix matrix::TKronT( const matrix& y ) M_CONST
  120. {
  121.    static char *mName = "kronT" ;
  122.    matFunc func( mName ) ; debugInfo( func ) ;
  123.    matrix z ;
  124.    z.kronOf( *this, y, TRANS12 ) ;
  125.    return z ;
  126. } // matrix::TKronT
  127.  
  128. matrix matrix::operator % ( const matrix& y ) M_CONST
  129. {
  130.    static char *mName = "op %" ;
  131.    matFunc func( mName ) ; debugInfo( func ) ;
  132.    matrix z ;
  133.    if ( isTrans() ) {
  134.       if ( y.isTrans() )
  135.          z.kronOf( *this, y, TRANS12 ) ;
  136.       else
  137.          z.kronOf( *this, y, TRANS1 ) ;
  138.    } else if ( y.isTrans() )
  139.       z.kronOf( *this, y, TRANS2 ) ;
  140.    else
  141.       z.kronOf( *this, y ) ;
  142.    return z ;
  143. } // matrix % (kron)
  144.