xref: /aosp_15_r20/external/eigen/test/permutationmatrices.cpp (revision bf2c37156dfe67e5dfebd6d394bad8b2ab5804d4)
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2009 Benoit Jacob <[email protected]>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #define TEST_ENABLE_TEMPORARY_TRACKING
11 
12 #include "main.h"
13 
14 using namespace std;
permutationmatrices(const MatrixType & m)15 template<typename MatrixType> void permutationmatrices(const MatrixType& m)
16 {
17   typedef typename MatrixType::Scalar Scalar;
18   enum { Rows = MatrixType::RowsAtCompileTime, Cols = MatrixType::ColsAtCompileTime,
19          Options = MatrixType::Options };
20   typedef PermutationMatrix<Rows> LeftPermutationType;
21   typedef Transpositions<Rows> LeftTranspositionsType;
22   typedef Matrix<int, Rows, 1> LeftPermutationVectorType;
23   typedef Map<LeftPermutationType> MapLeftPerm;
24   typedef PermutationMatrix<Cols> RightPermutationType;
25   typedef Transpositions<Cols> RightTranspositionsType;
26   typedef Matrix<int, Cols, 1> RightPermutationVectorType;
27   typedef Map<RightPermutationType> MapRightPerm;
28 
29   Index rows = m.rows();
30   Index cols = m.cols();
31 
32   MatrixType m_original = MatrixType::Random(rows,cols);
33   LeftPermutationVectorType lv;
34   randomPermutationVector(lv, rows);
35   LeftPermutationType lp(lv);
36   RightPermutationVectorType rv;
37   randomPermutationVector(rv, cols);
38   RightPermutationType rp(rv);
39   LeftTranspositionsType lt(lv);
40   RightTranspositionsType rt(rv);
41   MatrixType m_permuted = MatrixType::Random(rows,cols);
42 
43   VERIFY_EVALUATION_COUNT(m_permuted = lp * m_original * rp, 1); // 1 temp for sub expression "lp * m_original"
44 
45   for (int i=0; i<rows; i++)
46     for (int j=0; j<cols; j++)
47         VERIFY_IS_APPROX(m_permuted(lv(i),j), m_original(i,rv(j)));
48 
49   Matrix<Scalar,Rows,Rows> lm(lp);
50   Matrix<Scalar,Cols,Cols> rm(rp);
51 
52   VERIFY_IS_APPROX(m_permuted, lm*m_original*rm);
53 
54   m_permuted = m_original;
55   VERIFY_EVALUATION_COUNT(m_permuted = lp * m_permuted * rp, 1);
56   VERIFY_IS_APPROX(m_permuted, lm*m_original*rm);
57 
58   LeftPermutationType lpi;
59   lpi = lp.inverse();
60   VERIFY_IS_APPROX(lpi*m_permuted,lp.inverse()*m_permuted);
61 
62   VERIFY_IS_APPROX(lp.inverse()*m_permuted*rp.inverse(), m_original);
63   VERIFY_IS_APPROX(lv.asPermutation().inverse()*m_permuted*rv.asPermutation().inverse(), m_original);
64   VERIFY_IS_APPROX(MapLeftPerm(lv.data(),lv.size()).inverse()*m_permuted*MapRightPerm(rv.data(),rv.size()).inverse(), m_original);
65 
66   VERIFY((lp*lp.inverse()).toDenseMatrix().isIdentity());
67   VERIFY((lv.asPermutation()*lv.asPermutation().inverse()).toDenseMatrix().isIdentity());
68   VERIFY((MapLeftPerm(lv.data(),lv.size())*MapLeftPerm(lv.data(),lv.size()).inverse()).toDenseMatrix().isIdentity());
69 
70   LeftPermutationVectorType lv2;
71   randomPermutationVector(lv2, rows);
72   LeftPermutationType lp2(lv2);
73   Matrix<Scalar,Rows,Rows> lm2(lp2);
74   VERIFY_IS_APPROX((lp*lp2).toDenseMatrix().template cast<Scalar>(), lm*lm2);
75   VERIFY_IS_APPROX((lv.asPermutation()*lv2.asPermutation()).toDenseMatrix().template cast<Scalar>(), lm*lm2);
76   VERIFY_IS_APPROX((MapLeftPerm(lv.data(),lv.size())*MapLeftPerm(lv2.data(),lv2.size())).toDenseMatrix().template cast<Scalar>(), lm*lm2);
77 
78   LeftPermutationType identityp;
79   identityp.setIdentity(rows);
80   VERIFY_IS_APPROX(m_original, identityp*m_original);
81 
82   // check inplace permutations
83   m_permuted = m_original;
84   VERIFY_EVALUATION_COUNT(m_permuted.noalias()= lp.inverse() * m_permuted, 1); // 1 temp to allocate the mask
85   VERIFY_IS_APPROX(m_permuted, lp.inverse()*m_original);
86 
87   m_permuted = m_original;
88   VERIFY_EVALUATION_COUNT(m_permuted.noalias() = m_permuted * rp.inverse(), 1); // 1 temp to allocate the mask
89   VERIFY_IS_APPROX(m_permuted, m_original*rp.inverse());
90 
91   m_permuted = m_original;
92   VERIFY_EVALUATION_COUNT(m_permuted.noalias() = lp * m_permuted, 1); // 1 temp to allocate the mask
93   VERIFY_IS_APPROX(m_permuted, lp*m_original);
94 
95   m_permuted = m_original;
96   VERIFY_EVALUATION_COUNT(m_permuted.noalias() = m_permuted * rp, 1); // 1 temp to allocate the mask
97   VERIFY_IS_APPROX(m_permuted, m_original*rp);
98 
99   if(rows>1 && cols>1)
100   {
101     lp2 = lp;
102     Index i = internal::random<Index>(0, rows-1);
103     Index j;
104     do j = internal::random<Index>(0, rows-1); while(j==i);
105     lp2.applyTranspositionOnTheLeft(i, j);
106     lm = lp;
107     lm.row(i).swap(lm.row(j));
108     VERIFY_IS_APPROX(lm, lp2.toDenseMatrix().template cast<Scalar>());
109 
110     RightPermutationType rp2 = rp;
111     i = internal::random<Index>(0, cols-1);
112     do j = internal::random<Index>(0, cols-1); while(j==i);
113     rp2.applyTranspositionOnTheRight(i, j);
114     rm = rp;
115     rm.col(i).swap(rm.col(j));
116     VERIFY_IS_APPROX(rm, rp2.toDenseMatrix().template cast<Scalar>());
117   }
118 
119   {
120     // simple compilation check
121     Matrix<Scalar, Cols, Cols> A = rp;
122     Matrix<Scalar, Cols, Cols> B = rp.transpose();
123     VERIFY_IS_APPROX(A, B.transpose());
124   }
125 
126   m_permuted = m_original;
127   lp = lt;
128   rp = rt;
129   VERIFY_EVALUATION_COUNT(m_permuted = lt * m_permuted * rt, 1);
130   VERIFY_IS_APPROX(m_permuted, lp*m_original*rp.transpose());
131 
132   VERIFY_IS_APPROX(lt.inverse()*m_permuted*rt.inverse(), m_original);
133 
134   // Check inplace transpositions
135   m_permuted = m_original;
136   VERIFY_IS_APPROX(m_permuted = lt * m_permuted, lp * m_original);
137   m_permuted = m_original;
138   VERIFY_IS_APPROX(m_permuted = lt.inverse() * m_permuted, lp.inverse() * m_original);
139   m_permuted = m_original;
140   VERIFY_IS_APPROX(m_permuted = m_permuted * rt, m_original * rt);
141   m_permuted = m_original;
142   VERIFY_IS_APPROX(m_permuted = m_permuted * rt.inverse(), m_original * rt.inverse());
143 }
144 
145 template<typename T>
bug890()146 void bug890()
147 {
148   typedef Matrix<T, Dynamic, Dynamic> MatrixType;
149   typedef Matrix<T, Dynamic, 1> VectorType;
150   typedef Stride<Dynamic,Dynamic> S;
151   typedef Map<MatrixType, Aligned, S> MapType;
152   typedef PermutationMatrix<Dynamic> Perm;
153 
154   VectorType v1(2), v2(2), op(4), rhs(2);
155   v1 << 666,667;
156   op << 1,0,0,1;
157   rhs << 42,42;
158 
159   Perm P(2);
160   P.indices() << 1, 0;
161 
162   MapType(v1.data(),2,1,S(1,1)) = P * MapType(rhs.data(),2,1,S(1,1));
163   VERIFY_IS_APPROX(v1, (P * rhs).eval());
164 
165   MapType(v1.data(),2,1,S(1,1)) = P.inverse() * MapType(rhs.data(),2,1,S(1,1));
166   VERIFY_IS_APPROX(v1, (P.inverse() * rhs).eval());
167 }
168 
EIGEN_DECLARE_TEST(permutationmatrices)169 EIGEN_DECLARE_TEST(permutationmatrices)
170 {
171   for(int i = 0; i < g_repeat; i++) {
172     CALL_SUBTEST_1( permutationmatrices(Matrix<float, 1, 1>()) );
173     CALL_SUBTEST_2( permutationmatrices(Matrix3f()) );
174     CALL_SUBTEST_3( permutationmatrices(Matrix<double,3,3,RowMajor>()) );
175     CALL_SUBTEST_4( permutationmatrices(Matrix4d()) );
176     CALL_SUBTEST_5( permutationmatrices(Matrix<double,40,60>()) );
177     CALL_SUBTEST_6( permutationmatrices(Matrix<double,Dynamic,Dynamic,RowMajor>(internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
178     CALL_SUBTEST_7( permutationmatrices(MatrixXcf(internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
179   }
180   CALL_SUBTEST_5( bug890<double>() );
181 }
182