Re: [eigen] Fast QR for 2x2, 3x3 |
[ Thread Index |
Date Index
| More lists.tuxfamily.org/eigen Archives
]
Assuming m_q and m_r are 2x2 matrices to store the results, matrix is the
input 2x2 matrix which is to be factorized, this code gives the QR
factorisation:
if (matrix(1,0) == 0.0) // nothing to do
{
m_q = MatrixType::Identity();
m_r = matrix;
}
else if (matrix(0,0) == 0.0) // transposition of rows
{
m_q(0,0) = 0.0;
m_q(0,1) = 1.0;
m_q(1,0) = 1.0;
m_q(1,1) = 0.0;
m_r(0,0) = matrix(1,0);
m_r(0,1) = matrix(1,1);
m_r(1,0) = 0.0;
m_r(1,1) = matrix(0,1);
}
else
{
// avoid overflow. is this neccesary?
// if not, just use Scalar er = 1.0/matrix.block<2,1>(0,0).norm();
Scalar er;
if (matrix(0,0) >= matrix(1,0))
{
Scalar ba = matrix(1,0)/matrix(0,0);
er = 1.0/matrix(0,0) * sqrt(1+ba*ba);
}
else
{
Scalar ab = matrix(0,0)/matrix(1,0);
er = 1.0/matrix(1,0) * sqrt(1+ab*ab);
}
m_q(0,0) = er * matrix(0,0);
m_q(1,1) = er * matrix(0,0);
m_q(0,1) = -er * matrix(1,0);
m_q(1,0) = er * matrix(1,0);
m_r(0,0) = 1.0/er;
m_r(1,0) = 0.0;
m_r(0,1) = er * (matrix(0,0)*matrix(0,1) + matrix(1,0)*matrix(1,1));
m_r(1,1) = er * (matrix(0,0)*matrix(1,1) - matrix(1,0)*matrix(0,1));
}
A specialisation for 3x3 will follow in a few days.
-- Markus