home *** CD-ROM | disk | FTP | other *** search
/ io Programmo 21 / IOPROG_21.ISO / SOFT / LIBMAT.ZIP / SVDDEC.CPP < prev    next >
Encoding:
C/C++ Source or Header  |  1992-10-16  |  6.8 KB  |  284 lines
  1. /**************************************************/
  2. /*        svddec.c source for svdDec class        */
  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 "svddec.hpp"
  17.  
  18. /*******************************************************/
  19. /*                    svdDec methods                   */
  20. /*******************************************************/
  21.  
  22. void svdDec::setNames( void )
  23. {
  24.    S.name( "decSV" ) ;
  25.    V.name( "decV" ) ;
  26. } // svdDec setNames
  27.  
  28. #ifdef __cplusplus
  29. svdDec::svdDec( void ) : matDec()
  30. #else
  31. svdDec::svdDec( void ) : ()
  32. #endif
  33. {
  34.    setNames() ;
  35. } // svdDec( void )
  36.  
  37.  
  38. #ifdef __cplusplus
  39. svdDec::svdDec( const matrix& x ) : matDec(x)
  40. #else
  41. svdDec::svdDec( const matrix& x ) : (x)
  42. #endif
  43. {
  44.    INDEX nc = x.nCols() ;
  45.    S.reset( nc ) ;
  46.    V.reset( nc, nc ) ;
  47.    setNames() ;
  48. } // svdDec( matrix )
  49.  
  50. #ifdef __cplusplus
  51. svdDec::svdDec( const svdDec& sd ) : matDec( sd )
  52. #else
  53. svdDec::svdDec( const svdDec& sd ) : ( sd )
  54. #endif
  55. {
  56.    S = sd.S ;
  57.    V = sd.V ;
  58. } // svdDec( svdDec& )
  59.  
  60. void svdDec::operator = ( const svdDec& sd )
  61. {
  62.    static char *mName = "op = &" ;
  63.    matFunc func( mName ) ; debugInfo( func ) ;
  64.    matDec::operator = (sd) ;
  65.    S = sd.S ;
  66.    V = sd.V ;
  67. } // svdDec = svdDec
  68.  
  69. svdDec::~svdDec( void ) {} // ~svdDec
  70.  
  71. void svdDec::svd( matrix& u, matrix& sv, matrix& v )
  72. {
  73.    static char *mName = "svd" ;
  74.    matFunc func( mName ) ; debugInfo( func ) ;
  75.    if ( ok( mName ) ) {
  76.       u  = m ;
  77.       sv = S ;
  78.       v  = V ;
  79.    } //if
  80. } // svdDec svd
  81.  
  82. void svdDec::release( matrix& u, matrix& sv, matrix& v )
  83. {
  84.    static char *mName = "release" ;
  85.    matFunc func( mName ) ; debugInfo( func ) ;
  86.    if ( ok( mName ) ) {
  87.       u.capture( m ) ;
  88.       sv.capture( S ) ;
  89.       v.capture( V ) ;
  90.       clear() ;
  91.    } //if
  92. } // svdDec svd
  93.  
  94. void svdDec::clear( void )
  95. {
  96.    static char *mName = "clear" ;
  97.    matFunc func( mName ) ; debugInfo( func ) ;
  98.    matDec::clear() ;
  99.    S.clear() ;
  100.    V.clear() ;   
  101. } // svdDec clear
  102.  
  103. void svdDec::assign( const matrix& x )
  104. {
  105.    static char *mName = "assign" ;
  106.    matFunc func( mName ) ; debugInfo( func ) ;
  107.    matDec::assign( x ) ;
  108.    INDEX nc = x.nCols() ;
  109.    S.reset( nc ) ;
  110.    V.reset( nc, nc ) ;
  111. } // svdDec assign
  112.  
  113. void svdDec::capture( matrix& x )
  114. {
  115.    static char *mName = "assign" ;
  116.    matFunc func( mName ) ; debugInfo( func ) ;
  117.    matDec::capture( x ) ;
  118.    INDEX nc = m.nCols() ;
  119.    S.reset( nc ) ;
  120.    V.reset( nc, nc ) ;
  121. } // svdDec capture 
  122.  
  123. outFile& svdDec::info( outFile& f ) M_CONST
  124. {
  125.    return matDec::decInfo( f, "svdDec" ) ;
  126. } // qrhDec::info
  127.  
  128. inFile& svdDec::get( inFile& f )
  129. {
  130.    static char *mName = "get" ;
  131.    matFunc func( mName ) ; debugInfo( func ) ;
  132.    matDec::get(f) ;
  133.    INDEX n = m.nCols() ;
  134.    S.reset( n ) ;
  135.    V.reset( n, n ) ;
  136.    S.get(f) ;
  137.    V.get(f) ;
  138.    return f ;
  139. } // svdDec::get
  140.  
  141. outFile& svdDec::put( outFile& f ) M_CONST
  142. {
  143.    static char *mName = "put" ;
  144.    matFunc func( mName ) ; debugInfo( func ) ;
  145.    matDec::put(f) ;
  146.    S.put(f) ;
  147.    V.put(f) ;
  148.    return f ;
  149. } // svdDec::put
  150.  
  151. outFile& svdDec::print( char* decName, outFile& f ) M_CONST
  152. {
  153.    static char *mName = "print" ;
  154.    matFunc func( mName ) ; debugInfo( func ) ;
  155.    f.write( "SVD Decomposition : " ) ;
  156.    if ( decName != 0 )
  157.       f.write( decName ) ;
  158.    f.newLine(2) ;
  159.    matDec::print( "`U'" , f ) ;
  160.    S.print( "Singular Values : ", f ) ;
  161.    V.print( "`V' : ", f ) ;
  162.    return f ;
  163. } // svdDec::print
  164.  
  165. void svdDec::decompose( void )
  166. {
  167.    static char *mName = "decompose" ;
  168.    matFunc func( mName ) ; debugInfo( func ) ;
  169.    matError error ;
  170.    if ( !( status & ASSIGNED ) )
  171.       errorExit( mName, UNASS ) ;
  172.    if ( !( status & DECOMPOSED ) ) {
  173.       if ( !svdcmp( m, S, V, error ) )
  174.          errorExit( mName, error ) ;
  175.       status |= DECOMPOSED ;
  176.    } // if
  177. } // svdDec::decompose
  178.  
  179. void svdDec::zero( REAL min )
  180. {
  181.    static char *mName = "zero" ;
  182.    matFunc func( mName ) ; debugInfo( func ) ;
  183.    INDEX n = S.nRows(), i ;
  184.    if ( min <= 0.0 )
  185.       min = tol ;
  186.    if ( ok( mName ) ) {
  187.       for ( i = 1 ; i <= n ; i++ ) {
  188.          if ( S(i) < min )
  189.         S(i) = 0.0 ;
  190.       } // for
  191.    } // if
  192. } // svdDec:zero
  193.  
  194. void svdDec::sv( matrix& s )
  195. {
  196.    static char *mName = "sv" ;
  197.    matFunc func( mName ) ; debugInfo( func ) ;
  198.    if ( ok( mName ) ) 
  199.       s = S ;
  200. } // svdDec sv
  201.  
  202. void svdDec::setSV( matrix& s )
  203. {
  204.    static char *mName = "setSV" ;
  205.    matFunc func( mName ) ; debugInfo( func ) ;
  206.    if ( !( status & ASSIGNED ) )
  207.       errorExit( mName, UNASS ) ;
  208.    if ( S.nRows() != s.nRows() )
  209.       errorExit( mName, NEDIM ) ;
  210.    S = s ;
  211. } // svdDec sv
  212.  
  213. void svdDec::solve( matrix& b )
  214. {
  215.    static char *mName = "solve" ;
  216.    matFunc func( mName ) ; debugInfo( func ) ;
  217.    if ( ok( mName ) ) {
  218.       svdBackSub( m, S, V, b, b, error ) ;
  219.       if ( error )
  220.      errorExit( mName, error ) ;
  221.    } // else
  222. } // svdDec::solve
  223.  
  224. void svdDec::transSolve( matrix& b )
  225. {
  226.    static char *mName = "transSolve" ;
  227.    matFunc func( mName ) ; debugInfo( func ) ;
  228.    if ( ok( mName ) ) {
  229.       svdBackSub( V, S, m, b, b, error ) ;
  230.       if ( error )
  231.      errorExit( mName, error ) ;
  232.    } // else
  233. } // svdDec::transSolve
  234.  
  235. REAL svdDec::cond( void )
  236. {
  237.    static char *mName = "cond" ;
  238.    matFunc func( mName ) ; debugInfo( func ) ;
  239.    if ( ok( mName ) ) {
  240.       INDEX n = S.nRows() ;
  241.       REAL max = S(1) ;
  242.       REAL min = S(n) ;
  243.       if ( min == 0.0 )
  244.          condition = -1 ;
  245.       else
  246.          condition = max / min ;              
  247.       status |= CONDITION ;
  248.    } // if
  249.    return condition ;
  250. } // matDec::cond
  251.  
  252. INDEX svdDec::rank( void )
  253. {
  254.    static char *mName = "rank" ;
  255.    matFunc func( mName ) ; debugInfo( func ) ;
  256.    INDEX rankEst = 0 ;
  257.    if ( ok( mName ) ) {
  258.       for ( INDEX n = S.nRows() ; n ; n-- ) {
  259.          if ( S(n) > 0.0 )
  260.             rankEst++ ;
  261.       } // for
  262.    } // if
  263.    return rankEst ;
  264. } // svdDec rank
  265.  
  266. void svdDec::det( REAL& d1, REAL& d2 )
  267. {
  268.    static char *mName = "det" ;
  269.    matFunc func( mName ) ; debugInfo( func ) ;
  270.    if ( ok( mName ) ) {
  271.       if ( S( S.nRows() ) > 0 ) {
  272.          dProduct( S, det1, det2, 0.0, error ) ;
  273.          if ( error )
  274.             errorExit( mName, error ) ;
  275.       } else {
  276.          det1 = 0.0 ;
  277.          det2 = 0.0 ;
  278.       } // else
  279.       status |= DETERMINED ;
  280.    } // if
  281.    d1 = det1 ;
  282.    d2 = det2 ;   
  283. } // svdDec::det
  284.