ground/gcs/src/libs/eigen/test/vectorwiseop.cpp

   1 // This file is part of Eigen, a lightweight C++ template library
   2 // for linear algebra.
   3 //
   4 // Copyright (C) 2011 Benoit Jacob <jacob.benoit.1@gmail.com>
   5 // Copyright (C) 2015 Gael Guennebaud <gael.guennebaud@inria.fr>
   6 //
   7 // This Source Code Form is subject to the terms of the Mozilla
   8 // Public License v. 2.0. If a copy of the MPL was not distributed
   9 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
  10
  11 #define TEST_ENABLE_TEMPORARY_TRACKING
  12 #define EIGEN_NO_STATIC_ASSERT
  13
  14 #include "main.h"
  15
  16 template<typename ArrayType> void vectorwiseop_array(const ArrayType& m)
  17 {
  18   typedef typename ArrayType::Index Index;
  19   typedef typename ArrayType::Scalar Scalar;
  20   typedef Array<Scalar, ArrayType::RowsAtCompileTime, 1> ColVectorType;
  21   typedef Array<Scalar, 1, ArrayType::ColsAtCompileTime> RowVectorType;
  22
  23   Index rows = m.rows();
  24   Index cols = m.cols();
  25   Index r = internal::random<Index>(0, rows-1),
  26         c = internal::random<Index>(0, cols-1);
  27
  28   ArrayType m1 = ArrayType::Random(rows, cols),
  29             m2(rows, cols),
  30             m3(rows, cols);
  31
  32   ColVectorType colvec = ColVectorType::Random(rows);
  33   RowVectorType rowvec = RowVectorType::Random(cols);
  34
  35   // test addition
  36
  37   m2 = m1;
  38   m2.colwise() += colvec;
  39   VERIFY_IS_APPROX(m2, m1.colwise() + colvec);
  40   VERIFY_IS_APPROX(m2.col(c), m1.col(c) + colvec);
  41
  42   VERIFY_RAISES_ASSERT(m2.colwise() += colvec.transpose());
  43   VERIFY_RAISES_ASSERT(m1.colwise() + colvec.transpose());
  44
  45   m2 = m1;
  46   m2.rowwise() += rowvec;
  47   VERIFY_IS_APPROX(m2, m1.rowwise() + rowvec);
  48   VERIFY_IS_APPROX(m2.row(r), m1.row(r) + rowvec);
  49
  50   VERIFY_RAISES_ASSERT(m2.rowwise() += rowvec.transpose());
  51   VERIFY_RAISES_ASSERT(m1.rowwise() + rowvec.transpose());
  52
  53   // test substraction
  54
  55   m2 = m1;
  56   m2.colwise() -= colvec;
  57   VERIFY_IS_APPROX(m2, m1.colwise() - colvec);
  58   VERIFY_IS_APPROX(m2.col(c), m1.col(c) - colvec);
  59
  60   VERIFY_RAISES_ASSERT(m2.colwise() -= colvec.transpose());
  61   VERIFY_RAISES_ASSERT(m1.colwise() - colvec.transpose());
  62
  63   m2 = m1;
  64   m2.rowwise() -= rowvec;
  65   VERIFY_IS_APPROX(m2, m1.rowwise() - rowvec);
  66   VERIFY_IS_APPROX(m2.row(r), m1.row(r) - rowvec);
  67
  68   VERIFY_RAISES_ASSERT(m2.rowwise() -= rowvec.transpose());
  69   VERIFY_RAISES_ASSERT(m1.rowwise() - rowvec.transpose());
  70
  71   // test multiplication
  72
  73   m2 = m1;
  74   m2.colwise() *= colvec;
  75   VERIFY_IS_APPROX(m2, m1.colwise() * colvec);
  76   VERIFY_IS_APPROX(m2.col(c), m1.col(c) * colvec);
  77
  78   VERIFY_RAISES_ASSERT(m2.colwise() *= colvec.transpose());
  79   VERIFY_RAISES_ASSERT(m1.colwise() * colvec.transpose());
  80
  81   m2 = m1;
  82   m2.rowwise() *= rowvec;
  83   VERIFY_IS_APPROX(m2, m1.rowwise() * rowvec);
  84   VERIFY_IS_APPROX(m2.row(r), m1.row(r) * rowvec);
  85
  86   VERIFY_RAISES_ASSERT(m2.rowwise() *= rowvec.transpose());
  87   VERIFY_RAISES_ASSERT(m1.rowwise() * rowvec.transpose());
  88
  89   // test quotient
  90
  91   m2 = m1;
  92   m2.colwise() /= colvec;
  93   VERIFY_IS_APPROX(m2, m1.colwise() / colvec);
  94   VERIFY_IS_APPROX(m2.col(c), m1.col(c) / colvec);
  95
  96   VERIFY_RAISES_ASSERT(m2.colwise() /= colvec.transpose());
  97   VERIFY_RAISES_ASSERT(m1.colwise() / colvec.transpose());
  98
  99   m2 = m1;
 100   m2.rowwise() /= rowvec;
 101   VERIFY_IS_APPROX(m2, m1.rowwise() / rowvec);
 102   VERIFY_IS_APPROX(m2.row(r), m1.row(r) / rowvec);
 103
 104   VERIFY_RAISES_ASSERT(m2.rowwise() /= rowvec.transpose());
 105   VERIFY_RAISES_ASSERT(m1.rowwise() / rowvec.transpose());
 106
 107   m2 = m1;
 108   // yes, there might be an aliasing issue there but ".rowwise() /="
 109   // is supposed to evaluate " m2.colwise().sum()" into a temporary to avoid
 110   // evaluating the reduction multiple times
 111   if(ArrayType::RowsAtCompileTime>2 || ArrayType::RowsAtCompileTime==Dynamic)
 112   {
 113     m2.rowwise() /= m2.colwise().sum();
 114     VERIFY_IS_APPROX(m2, m1.rowwise() / m1.colwise().sum());
 115   }
 116
 117   // all/any
 118   Array<bool,Dynamic,Dynamic> mb(rows,cols);
 119   mb = (m1.real()<=0.7).colwise().all();
 120   VERIFY( (mb.col(c) == (m1.real().col(c)<=0.7).all()).all() );
 121   mb = (m1.real()<=0.7).rowwise().all();
 122   VERIFY( (mb.row(r) == (m1.real().row(r)<=0.7).all()).all() );
 123
 124   mb = (m1.real()>=0.7).colwise().any();
 125   VERIFY( (mb.col(c) == (m1.real().col(c)>=0.7).any()).all() );
 126   mb = (m1.real()>=0.7).rowwise().any();
 127   VERIFY( (mb.row(r) == (m1.real().row(r)>=0.7).any()).all() );
 128 }
 129
 130 template<typename MatrixType> void vectorwiseop_matrix(const MatrixType& m)
 131 {
 132   typedef typename MatrixType::Index Index;
 133   typedef typename MatrixType::Scalar Scalar;
 134   typedef typename NumTraits<Scalar>::Real RealScalar;
 135   typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> ColVectorType;
 136   typedef Matrix<Scalar, 1, MatrixType::ColsAtCompileTime> RowVectorType;
 137   typedef Matrix<RealScalar, MatrixType::RowsAtCompileTime, 1> RealColVectorType;
 138   typedef Matrix<RealScalar, 1, MatrixType::ColsAtCompileTime> RealRowVectorType;
 139
 140   Index rows = m.rows();
 141   Index cols = m.cols();
 142   Index r = internal::random<Index>(0, rows-1),
 143         c = internal::random<Index>(0, cols-1);
 144
 145   MatrixType m1 = MatrixType::Random(rows, cols),
 146             m2(rows, cols),
 147             m3(rows, cols);
 148
 149   ColVectorType colvec = ColVectorType::Random(rows);
 150   RowVectorType rowvec = RowVectorType::Random(cols);
 151   RealColVectorType rcres;
 152   RealRowVectorType rrres;
 153
 154   // test addition
 155
 156   m2 = m1;
 157   m2.colwise() += colvec;
 158   VERIFY_IS_APPROX(m2, m1.colwise() + colvec);
 159   VERIFY_IS_APPROX(m2.col(c), m1.col(c) + colvec);
 160
 161   if(rows>1)
 162   {
 163     VERIFY_RAISES_ASSERT(m2.colwise() += colvec.transpose());
 164     VERIFY_RAISES_ASSERT(m1.colwise() + colvec.transpose());
 165   }
 166
 167   m2 = m1;
 168   m2.rowwise() += rowvec;
 169   VERIFY_IS_APPROX(m2, m1.rowwise() + rowvec);
 170   VERIFY_IS_APPROX(m2.row(r), m1.row(r) + rowvec);
 171
 172   if(cols>1)
 173   {
 174     VERIFY_RAISES_ASSERT(m2.rowwise() += rowvec.transpose());
 175     VERIFY_RAISES_ASSERT(m1.rowwise() + rowvec.transpose());
 176   }
 177
 178   // test substraction
 179
 180   m2 = m1;
 181   m2.colwise() -= colvec;
 182   VERIFY_IS_APPROX(m2, m1.colwise() - colvec);
 183   VERIFY_IS_APPROX(m2.col(c), m1.col(c) - colvec);
 184
 185   if(rows>1)
 186   {
 187     VERIFY_RAISES_ASSERT(m2.colwise() -= colvec.transpose());
 188     VERIFY_RAISES_ASSERT(m1.colwise() - colvec.transpose());
 189   }
 190
 191   m2 = m1;
 192   m2.rowwise() -= rowvec;
 193   VERIFY_IS_APPROX(m2, m1.rowwise() - rowvec);
 194   VERIFY_IS_APPROX(m2.row(r), m1.row(r) - rowvec);
 195
 196   if(cols>1)
 197   {
 198     VERIFY_RAISES_ASSERT(m2.rowwise() -= rowvec.transpose());
 199     VERIFY_RAISES_ASSERT(m1.rowwise() - rowvec.transpose());
 200   }
 201
 202   // test norm
 203   rrres = m1.colwise().norm();
 204   VERIFY_IS_APPROX(rrres(c), m1.col(c).norm());
 205   rcres = m1.rowwise().norm();
 206   VERIFY_IS_APPROX(rcres(r), m1.row(r).norm());
 207
 208   VERIFY_IS_APPROX(m1.cwiseAbs().colwise().sum(), m1.colwise().template lpNorm<1>());
 209   VERIFY_IS_APPROX(m1.cwiseAbs().rowwise().sum(), m1.rowwise().template lpNorm<1>());
 210   VERIFY_IS_APPROX(m1.cwiseAbs().colwise().maxCoeff(), m1.colwise().template lpNorm<Infinity>());
 211   VERIFY_IS_APPROX(m1.cwiseAbs().rowwise().maxCoeff(), m1.rowwise().template lpNorm<Infinity>());
 212
 213   // regression for bug 1158
 214   VERIFY_IS_APPROX(m1.cwiseAbs().colwise().sum().x(), m1.col(0).cwiseAbs().sum());
 215
 216   // test normalized
 217   m2 = m1.colwise().normalized();
 218   VERIFY_IS_APPROX(m2.col(c), m1.col(c).normalized());
 219   m2 = m1.rowwise().normalized();
 220   VERIFY_IS_APPROX(m2.row(r), m1.row(r).normalized());
 221
 222   // test normalize
 223   m2 = m1;
 224   m2.colwise().normalize();
 225   VERIFY_IS_APPROX(m2.col(c), m1.col(c).normalized());
 226   m2 = m1;
 227   m2.rowwise().normalize();
 228   VERIFY_IS_APPROX(m2.row(r), m1.row(r).normalized());
 229
 230   // test with partial reduction of products
 231   Matrix<Scalar,MatrixType::RowsAtCompileTime,MatrixType::RowsAtCompileTime> m1m1 = m1 * m1.transpose();
 232   VERIFY_IS_APPROX( (m1 * m1.transpose()).colwise().sum(), m1m1.colwise().sum());
 233   Matrix<Scalar,1,MatrixType::RowsAtCompileTime> tmp(rows);
 234   VERIFY_EVALUATION_COUNT( tmp = (m1 * m1.transpose()).colwise().sum(), 1);
 235
 236   m2 = m1.rowwise() - (m1.colwise().sum()/RealScalar(m1.rows())).eval();
 237   m1 = m1.rowwise() - (m1.colwise().sum()/RealScalar(m1.rows()));
 238   VERIFY_IS_APPROX( m1, m2 );
 239   VERIFY_EVALUATION_COUNT( m2 = (m1.rowwise() - m1.colwise().sum()/RealScalar(m1.rows())), (MatrixType::RowsAtCompileTime!=1 ? 1 : 0) );
 240 }
 241
 242 void test_vectorwiseop()
 243 {
 244   CALL_SUBTEST_1( vectorwiseop_array(Array22cd()) );
 245   CALL_SUBTEST_2( vectorwiseop_array(Array<double, 3, 2>()) );
 246   CALL_SUBTEST_3( vectorwiseop_array(ArrayXXf(3, 4)) );
 247   CALL_SUBTEST_4( vectorwiseop_matrix(Matrix4cf()) );
 248   CALL_SUBTEST_5( vectorwiseop_matrix(Matrix<float,4,5>()) );
 249   CALL_SUBTEST_6( vectorwiseop_matrix(MatrixXd(internal::random<int>(1,EIGEN_TEST_MAX_SIZE), internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
 250   CALL_SUBTEST_7( vectorwiseop_matrix(VectorXd(internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
 251   CALL_SUBTEST_7( vectorwiseop_matrix(RowVectorXd(internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
 252 }