1 MatrixXf A
= MatrixXf::Random(4,4);
2 MatrixXf B
= MatrixXf::Random(4,4);
3 RealQZ
<MatrixXf
> qz(4); // preallocate space for 4x4 matrices
4 qz
.compute(A
,B
); // A = Q S Z, B = Q T Z
6 // print original matrices and result of decomposition
7 cout
<< "A:\n" << A
<< "\n" << "B:\n" << B
<< "\n";
8 cout
<< "S:\n" << qz
.matrixS() << "\n" << "T:\n" << qz
.matrixT() << "\n";
9 cout
<< "Q:\n" << qz
.matrixQ() << "\n" << "Z:\n" << qz
.matrixZ() << "\n";
13 << "\n|A-QSZ|: " << (A
-qz
.matrixQ()*qz
.matrixS()*qz
.matrixZ()).norm()
14 << ", |B-QTZ|: " << (B
-qz
.matrixQ()*qz
.matrixT()*qz
.matrixZ()).norm()
15 << "\n|QQ* - I|: " << (qz
.matrixQ()*qz
.matrixQ().adjoint() - MatrixXf::Identity(4,4)).norm()
16 << ", |ZZ* - I|: " << (qz
.matrixZ()*qz
.matrixZ().adjoint() - MatrixXf::Identity(4,4)).norm()