[eigen] pseudo-inverse

[ Thread Index | Date Index | More lists.tuxfamily.org/eigen Archives ]


hello -- here is some code i added to the SVD.h and FullPivHouseholderQR.h 

not sure if it's useful or if the qr is correct ... the svd seems to work. 

b.a.


void pinv( MatrixType& pinvmat)  
    {
      ei_assert(m_isInitialized && "SVD is not initialized.");
      double  pinvtoler=1..e-6;
      SingularValuesType m_sigma_inv=m_sigma;
  
      for ( long i=0; i<m_workMatrix.cols(); ++i)
{   
 if ( m_sigma(i) > pinvtoler )  // FIXME -- check tolerances against matlab pinv
   m_sigma_inv(i)=1.0/m_sigma(i);
 else m_sigma_inv(i)=0;
}
      pinvmat= (m_matV*m_sigma_inv.asDiagonal()*m_matU.transpose()); 
    }
    

    void pinv(MatrixType& pinvmat)  
    {
      ei_assert(m_isInitialized && "FullPivHouseholderQR is not initialized.");
     // compute R matrix 
      MatrixType temp=matrixQ().transpose()*m_qr;
     // use the fact that A*A = R*R  to get  (A*A)_inv * A*
      temp=(temp.transpose()*temp).inverse();
      pinvmat=temp*m_qr.transpose();
    }






--
ß®∫∆π 



Mail converted by MHonArc 2.6.19+ http://listengine.tuxfamily.org/