cart-elc

Source code for CART-ELC
git clone git://git.laack.co/cart-elc.git
Log | Files | Refs | README | LICENSE

product_small.cpp (12656B)


      1 // This file is part of Eigen, a lightweight C++ template library
      2 // for linear algebra.
      3 //
      4 // Copyright (C) 2006-2008 Benoit Jacob <jacob.benoit.1@gmail.com>
      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 EIGEN_NO_STATIC_ASSERT
     11 #include "product.h"
     12 #include <Eigen/LU>
     13 
     14 // regression test for bug 447
     15 template<int>
     16 void product1x1()
     17 {
     18   Matrix<float,1,3> matAstatic;
     19   Matrix<float,3,1> matBstatic;
     20   matAstatic.setRandom();
     21   matBstatic.setRandom();
     22   VERIFY_IS_APPROX( (matAstatic * matBstatic).coeff(0,0), 
     23                     matAstatic.cwiseProduct(matBstatic.transpose()).sum() );
     24 
     25   MatrixXf matAdynamic(1,3);
     26   MatrixXf matBdynamic(3,1);
     27   matAdynamic.setRandom();
     28   matBdynamic.setRandom();
     29   VERIFY_IS_APPROX( (matAdynamic * matBdynamic).coeff(0,0), 
     30                     matAdynamic.cwiseProduct(matBdynamic.transpose()).sum() );
     31 }
     32 
     33 template<typename TC, typename TA, typename TB>
     34 const TC& ref_prod(TC &C, const TA &A, const TB &B)
     35 {
     36   for(Index i=0;i<C.rows();++i)
     37     for(Index j=0;j<C.cols();++j)
     38       for(Index k=0;k<A.cols();++k)
     39         C.coeffRef(i,j) += A.coeff(i,k) * B.coeff(k,j);
     40   return C;
     41 }
     42 
     43 template<typename T, int Rows, int Cols, int Depth, int OC, int OA, int OB>
     44 typename internal::enable_if<! ( (Rows ==1&&Depth!=1&&OA==ColMajor)
     45                               || (Depth==1&&Rows !=1&&OA==RowMajor)
     46                               || (Cols ==1&&Depth!=1&&OB==RowMajor)
     47                               || (Depth==1&&Cols !=1&&OB==ColMajor)
     48                               || (Rows ==1&&Cols !=1&&OC==ColMajor)
     49                               || (Cols ==1&&Rows !=1&&OC==RowMajor)),void>::type
     50 test_lazy_single(int rows, int cols, int depth)
     51 {
     52   Matrix<T,Rows,Depth,OA> A(rows,depth); A.setRandom();
     53   Matrix<T,Depth,Cols,OB> B(depth,cols); B.setRandom();
     54   Matrix<T,Rows,Cols,OC>  C(rows,cols);  C.setRandom();
     55   Matrix<T,Rows,Cols,OC>  D(C);
     56   VERIFY_IS_APPROX(C+=A.lazyProduct(B), ref_prod(D,A,B));
     57 }
     58 
     59 void test_dynamic_bool()
     60 {
     61   int rows = internal::random<int>(1,64);
     62   int cols = internal::random<int>(1,64);
     63   int depth = internal::random<int>(1,65);
     64 
     65   typedef Matrix<bool,Dynamic,Dynamic> MatrixX;
     66   MatrixX A(rows,depth); A.setRandom();
     67   MatrixX B(depth,cols); B.setRandom();
     68   MatrixX C(rows,cols);  C.setRandom();
     69   MatrixX D(C);
     70   for(Index i=0;i<C.rows();++i)
     71     for(Index j=0;j<C.cols();++j)
     72       for(Index k=0;k<A.cols();++k)
     73        D.coeffRef(i,j) |= A.coeff(i,k) & B.coeff(k,j);
     74   C += A * B;
     75   VERIFY_IS_EQUAL(C, D);
     76 
     77   MatrixX E = B.transpose();
     78   for(Index i=0;i<B.rows();++i)
     79     for(Index j=0;j<B.cols();++j)
     80       VERIFY_IS_EQUAL(B(i,j), E(j,i));
     81 }
     82 
     83 template<typename T, int Rows, int Cols, int Depth, int OC, int OA, int OB>
     84 typename internal::enable_if<  ( (Rows ==1&&Depth!=1&&OA==ColMajor)
     85                               || (Depth==1&&Rows !=1&&OA==RowMajor)
     86                               || (Cols ==1&&Depth!=1&&OB==RowMajor)
     87                               || (Depth==1&&Cols !=1&&OB==ColMajor)
     88                               || (Rows ==1&&Cols !=1&&OC==ColMajor)
     89                               || (Cols ==1&&Rows !=1&&OC==RowMajor)),void>::type
     90 test_lazy_single(int, int, int)
     91 {
     92 }
     93 
     94 template<typename T, int Rows, int Cols, int Depth>
     95 void test_lazy_all_layout(int rows=Rows, int cols=Cols, int depth=Depth)
     96 {
     97   CALL_SUBTEST(( test_lazy_single<T,Rows,Cols,Depth,ColMajor,ColMajor,ColMajor>(rows,cols,depth) ));
     98   CALL_SUBTEST(( test_lazy_single<T,Rows,Cols,Depth,RowMajor,ColMajor,ColMajor>(rows,cols,depth) ));
     99   CALL_SUBTEST(( test_lazy_single<T,Rows,Cols,Depth,ColMajor,RowMajor,ColMajor>(rows,cols,depth) ));
    100   CALL_SUBTEST(( test_lazy_single<T,Rows,Cols,Depth,RowMajor,RowMajor,ColMajor>(rows,cols,depth) ));
    101   CALL_SUBTEST(( test_lazy_single<T,Rows,Cols,Depth,ColMajor,ColMajor,RowMajor>(rows,cols,depth) ));
    102   CALL_SUBTEST(( test_lazy_single<T,Rows,Cols,Depth,RowMajor,ColMajor,RowMajor>(rows,cols,depth) ));
    103   CALL_SUBTEST(( test_lazy_single<T,Rows,Cols,Depth,ColMajor,RowMajor,RowMajor>(rows,cols,depth) ));
    104   CALL_SUBTEST(( test_lazy_single<T,Rows,Cols,Depth,RowMajor,RowMajor,RowMajor>(rows,cols,depth) ));
    105 }  
    106 
    107 template<typename T>
    108 void test_lazy_l1()
    109 {
    110   int rows = internal::random<int>(1,12);
    111   int cols = internal::random<int>(1,12);
    112   int depth = internal::random<int>(1,12);
    113 
    114   // Inner
    115   CALL_SUBTEST(( test_lazy_all_layout<T,1,1,1>() ));
    116   CALL_SUBTEST(( test_lazy_all_layout<T,1,1,2>() ));
    117   CALL_SUBTEST(( test_lazy_all_layout<T,1,1,3>() ));
    118   CALL_SUBTEST(( test_lazy_all_layout<T,1,1,8>() ));
    119   CALL_SUBTEST(( test_lazy_all_layout<T,1,1,9>() ));
    120   CALL_SUBTEST(( test_lazy_all_layout<T,1,1,-1>(1,1,depth) ));
    121 
    122   // Outer
    123   CALL_SUBTEST(( test_lazy_all_layout<T,2,1,1>() ));
    124   CALL_SUBTEST(( test_lazy_all_layout<T,1,2,1>() ));
    125   CALL_SUBTEST(( test_lazy_all_layout<T,2,2,1>() ));
    126   CALL_SUBTEST(( test_lazy_all_layout<T,3,3,1>() ));
    127   CALL_SUBTEST(( test_lazy_all_layout<T,4,4,1>() ));
    128   CALL_SUBTEST(( test_lazy_all_layout<T,4,8,1>() ));
    129   CALL_SUBTEST(( test_lazy_all_layout<T,4,-1,1>(4,cols) ));
    130   CALL_SUBTEST(( test_lazy_all_layout<T,7,-1,1>(7,cols) ));
    131   CALL_SUBTEST(( test_lazy_all_layout<T,-1,8,1>(rows) ));
    132   CALL_SUBTEST(( test_lazy_all_layout<T,-1,3,1>(rows) ));
    133   CALL_SUBTEST(( test_lazy_all_layout<T,-1,-1,1>(rows,cols) ));
    134 }
    135 
    136 template<typename T>
    137 void test_lazy_l2()
    138 {
    139   int rows = internal::random<int>(1,12);
    140   int cols = internal::random<int>(1,12);
    141   int depth = internal::random<int>(1,12);
    142 
    143   // mat-vec
    144   CALL_SUBTEST(( test_lazy_all_layout<T,2,1,2>() ));
    145   CALL_SUBTEST(( test_lazy_all_layout<T,2,1,4>() ));
    146   CALL_SUBTEST(( test_lazy_all_layout<T,4,1,2>() ));
    147   CALL_SUBTEST(( test_lazy_all_layout<T,4,1,4>() ));
    148   CALL_SUBTEST(( test_lazy_all_layout<T,5,1,4>() ));
    149   CALL_SUBTEST(( test_lazy_all_layout<T,4,1,5>() ));
    150   CALL_SUBTEST(( test_lazy_all_layout<T,4,1,6>() ));
    151   CALL_SUBTEST(( test_lazy_all_layout<T,6,1,4>() ));
    152   CALL_SUBTEST(( test_lazy_all_layout<T,8,1,8>() ));
    153   CALL_SUBTEST(( test_lazy_all_layout<T,-1,1,4>(rows) ));
    154   CALL_SUBTEST(( test_lazy_all_layout<T,4,1,-1>(4,1,depth) ));
    155   CALL_SUBTEST(( test_lazy_all_layout<T,-1,1,-1>(rows,1,depth) ));
    156 
    157   // vec-mat
    158   CALL_SUBTEST(( test_lazy_all_layout<T,1,2,2>() ));
    159   CALL_SUBTEST(( test_lazy_all_layout<T,1,2,4>() ));
    160   CALL_SUBTEST(( test_lazy_all_layout<T,1,4,2>() ));
    161   CALL_SUBTEST(( test_lazy_all_layout<T,1,4,4>() ));
    162   CALL_SUBTEST(( test_lazy_all_layout<T,1,5,4>() ));
    163   CALL_SUBTEST(( test_lazy_all_layout<T,1,4,5>() ));
    164   CALL_SUBTEST(( test_lazy_all_layout<T,1,4,6>() ));
    165   CALL_SUBTEST(( test_lazy_all_layout<T,1,6,4>() ));
    166   CALL_SUBTEST(( test_lazy_all_layout<T,1,8,8>() ));
    167   CALL_SUBTEST(( test_lazy_all_layout<T,1,-1, 4>(1,cols) ));
    168   CALL_SUBTEST(( test_lazy_all_layout<T,1, 4,-1>(1,4,depth) ));
    169   CALL_SUBTEST(( test_lazy_all_layout<T,1,-1,-1>(1,cols,depth) ));
    170 }
    171 
    172 template<typename T>
    173 void test_lazy_l3()
    174 {
    175   int rows = internal::random<int>(1,12);
    176   int cols = internal::random<int>(1,12);
    177   int depth = internal::random<int>(1,12);
    178   // mat-mat
    179   CALL_SUBTEST(( test_lazy_all_layout<T,2,4,2>() ));
    180   CALL_SUBTEST(( test_lazy_all_layout<T,2,6,4>() ));
    181   CALL_SUBTEST(( test_lazy_all_layout<T,4,3,2>() ));
    182   CALL_SUBTEST(( test_lazy_all_layout<T,4,8,4>() ));
    183   CALL_SUBTEST(( test_lazy_all_layout<T,5,6,4>() ));
    184   CALL_SUBTEST(( test_lazy_all_layout<T,4,2,5>() ));
    185   CALL_SUBTEST(( test_lazy_all_layout<T,4,7,6>() ));
    186   CALL_SUBTEST(( test_lazy_all_layout<T,6,8,4>() ));
    187   CALL_SUBTEST(( test_lazy_all_layout<T,8,3,8>() ));
    188   CALL_SUBTEST(( test_lazy_all_layout<T,-1,6,4>(rows) ));
    189   CALL_SUBTEST(( test_lazy_all_layout<T,4,3,-1>(4,3,depth) ));
    190   CALL_SUBTEST(( test_lazy_all_layout<T,-1,6,-1>(rows,6,depth) ));
    191   CALL_SUBTEST(( test_lazy_all_layout<T,8,2,2>() ));
    192   CALL_SUBTEST(( test_lazy_all_layout<T,5,2,4>() ));
    193   CALL_SUBTEST(( test_lazy_all_layout<T,4,4,2>() ));
    194   CALL_SUBTEST(( test_lazy_all_layout<T,8,4,4>() ));
    195   CALL_SUBTEST(( test_lazy_all_layout<T,6,5,4>() ));
    196   CALL_SUBTEST(( test_lazy_all_layout<T,4,4,5>() ));
    197   CALL_SUBTEST(( test_lazy_all_layout<T,3,4,6>() ));
    198   CALL_SUBTEST(( test_lazy_all_layout<T,2,6,4>() ));
    199   CALL_SUBTEST(( test_lazy_all_layout<T,7,8,8>() ));
    200   CALL_SUBTEST(( test_lazy_all_layout<T,8,-1, 4>(8,cols) ));
    201   CALL_SUBTEST(( test_lazy_all_layout<T,3, 4,-1>(3,4,depth) ));
    202   CALL_SUBTEST(( test_lazy_all_layout<T,4,-1,-1>(4,cols,depth) ));
    203 }
    204 
    205 template<typename T,int N,int M,int K>
    206 void test_linear_but_not_vectorizable()
    207 {
    208   // Check tricky cases for which the result of the product is a vector and thus must exhibit the LinearBit flag,
    209   // but is not vectorizable along the linear dimension.
    210   Index n = N==Dynamic ? internal::random<Index>(1,32) : N;
    211   Index m = M==Dynamic ? internal::random<Index>(1,32) : M;
    212   Index k = K==Dynamic ? internal::random<Index>(1,32) : K;
    213 
    214   {
    215     Matrix<T,N,M+1> A; A.setRandom(n,m+1);
    216     Matrix<T,M*2,K> B; B.setRandom(m*2,k);
    217     Matrix<T,1,K> C;
    218     Matrix<T,1,K> R;
    219 
    220     C.noalias() = A.template topLeftCorner<1,M>() * (B.template topRows<M>()+B.template bottomRows<M>());
    221     R.noalias() = A.template topLeftCorner<1,M>() * (B.template topRows<M>()+B.template bottomRows<M>()).eval();
    222     VERIFY_IS_APPROX(C,R);
    223   }
    224 
    225   {
    226     Matrix<T,M+1,N,RowMajor> A; A.setRandom(m+1,n);
    227     Matrix<T,K,M*2,RowMajor> B; B.setRandom(k,m*2);
    228     Matrix<T,K,1> C;
    229     Matrix<T,K,1> R;
    230 
    231     C.noalias() = (B.template leftCols<M>()+B.template rightCols<M>())        * A.template topLeftCorner<M,1>();
    232     R.noalias() = (B.template leftCols<M>()+B.template rightCols<M>()).eval() * A.template topLeftCorner<M,1>();
    233     VERIFY_IS_APPROX(C,R);
    234   }
    235 }
    236 
    237 template<int Rows>
    238 void bug_1311()
    239 {
    240   Matrix< double, Rows, 2 > A;  A.setRandom();
    241   Vector2d b = Vector2d::Random() ;
    242   Matrix<double,Rows,1> res;
    243   res.noalias() = 1. * (A * b);
    244   VERIFY_IS_APPROX(res, A*b);
    245   res.noalias() = 1.*A * b;
    246   VERIFY_IS_APPROX(res, A*b);
    247   res.noalias() = (1.*A).lazyProduct(b);
    248   VERIFY_IS_APPROX(res, A*b);
    249   res.noalias() = (1.*A).lazyProduct(1.*b);
    250   VERIFY_IS_APPROX(res, A*b);
    251   res.noalias() = (A).lazyProduct(1.*b);
    252   VERIFY_IS_APPROX(res, A*b);
    253 }
    254 
    255 template<int>
    256 void product_small_regressions()
    257 {
    258   {
    259     // test compilation of (outer_product) * vector
    260     Vector3f v = Vector3f::Random();
    261     VERIFY_IS_APPROX( (v * v.transpose()) * v, (v * v.transpose()).eval() * v);
    262   }
    263   
    264   {
    265     // regression test for pull-request #93
    266     Eigen::Matrix<double, 1, 1> A;  A.setRandom();
    267     Eigen::Matrix<double, 18, 1> B; B.setRandom();
    268     Eigen::Matrix<double, 1, 18> C; C.setRandom();
    269     VERIFY_IS_APPROX(B * A.inverse(), B * A.inverse()[0]);
    270     VERIFY_IS_APPROX(A.inverse() * C, A.inverse()[0] * C);
    271   }
    272 
    273   {
    274     Eigen::Matrix<double, 10, 10> A, B, C;
    275     A.setRandom();
    276     C = A;
    277     for(int k=0; k<79; ++k)
    278       C = C * A;
    279     B.noalias() = (((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A)) * ((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A)))
    280                 * (((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A)) * ((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A))*((A*A)*(A*A)));
    281     VERIFY_IS_APPROX(B,C);
    282   }
    283 }
    284 
    285 EIGEN_DECLARE_TEST(product_small)
    286 {
    287   for(int i = 0; i < g_repeat; i++) {
    288     CALL_SUBTEST_1( product(Matrix<float, 3, 2>()) );
    289     CALL_SUBTEST_2( product(Matrix<int, 3, 17>()) );
    290     CALL_SUBTEST_8( product(Matrix<double, 3, 17>()) );
    291     CALL_SUBTEST_3( product(Matrix3d()) );
    292     CALL_SUBTEST_4( product(Matrix4d()) );
    293     CALL_SUBTEST_5( product(Matrix4f()) );
    294     CALL_SUBTEST_6( product1x1<0>() );
    295 
    296     CALL_SUBTEST_11( test_lazy_l1<float>() );
    297     CALL_SUBTEST_12( test_lazy_l2<float>() );
    298     CALL_SUBTEST_13( test_lazy_l3<float>() );
    299 
    300     CALL_SUBTEST_21( test_lazy_l1<double>() );
    301     CALL_SUBTEST_22( test_lazy_l2<double>() );
    302     CALL_SUBTEST_23( test_lazy_l3<double>() );
    303 
    304     CALL_SUBTEST_31( test_lazy_l1<std::complex<float> >() );
    305     CALL_SUBTEST_32( test_lazy_l2<std::complex<float> >() );
    306     CALL_SUBTEST_33( test_lazy_l3<std::complex<float> >() );
    307 
    308     CALL_SUBTEST_41( test_lazy_l1<std::complex<double> >() );
    309     CALL_SUBTEST_42( test_lazy_l2<std::complex<double> >() );
    310     CALL_SUBTEST_43( test_lazy_l3<std::complex<double> >() );
    311 
    312     CALL_SUBTEST_7(( test_linear_but_not_vectorizable<float,2,1,Dynamic>() ));
    313     CALL_SUBTEST_7(( test_linear_but_not_vectorizable<float,3,1,Dynamic>() ));
    314     CALL_SUBTEST_7(( test_linear_but_not_vectorizable<float,2,1,16>() ));
    315 
    316     CALL_SUBTEST_6( bug_1311<3>() );
    317     CALL_SUBTEST_6( bug_1311<5>() );
    318 
    319     CALL_SUBTEST_9( test_dynamic_bool() );
    320   }
    321 
    322   CALL_SUBTEST_6( product_small_regressions<0>() );
    323 }