array_for_matrix.cpp
Go to the documentation of this file.
00001 // This file is part of Eigen, a lightweight C++ template library
00002 // for linear algebra.
00003 //
00004 // Copyright (C) 2008-2009 Gael Guennebaud <gael.guennebaud@inria.fr>
00005 //
00006 // Eigen is free software; you can redistribute it and/or
00007 // modify it under the terms of the GNU Lesser General Public
00008 // License as published by the Free Software Foundation; either
00009 // version 3 of the License, or (at your option) any later version.
00010 //
00011 // Alternatively, you can redistribute it and/or
00012 // modify it under the terms of the GNU General Public License as
00013 // published by the Free Software Foundation; either version 2 of
00014 // the License, or (at your option) any later version.
00015 //
00016 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
00017 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
00018 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
00019 // GNU General Public License for more details.
00020 //
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License and a copy of the GNU General Public License along with
00023 // Eigen. If not, see <http://www.gnu.org/licenses/>.
00024 
00025 #include "main.h"
00026 
00027 template<typename MatrixType> void array_for_matrix(const MatrixType& m)
00028 {
00029   typedef typename MatrixType::Index Index;
00030   typedef typename MatrixType::Scalar Scalar;
00031   typedef typename NumTraits<Scalar>::Real RealScalar;
00032   typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> ColVectorType;
00033   typedef Matrix<Scalar, 1, MatrixType::ColsAtCompileTime> RowVectorType; 
00034 
00035   Index rows = m.rows();
00036   Index cols = m.cols();
00037 
00038   MatrixType m1 = MatrixType::Random(rows, cols),
00039              m2 = MatrixType::Random(rows, cols),
00040              m3(rows, cols);
00041 
00042   ColVectorType cv1 = ColVectorType::Random(rows);
00043   RowVectorType rv1 = RowVectorType::Random(cols);
00044 
00045   Scalar  s1 = internal::random<Scalar>(),
00046           s2 = internal::random<Scalar>();
00047 
00048   // scalar addition
00049   VERIFY_IS_APPROX(m1.array() + s1, s1 + m1.array());
00050   VERIFY_IS_APPROX((m1.array() + s1).matrix(), MatrixType::Constant(rows,cols,s1) + m1);
00051   VERIFY_IS_APPROX(((m1*Scalar(2)).array() - s2).matrix(), (m1+m1) - MatrixType::Constant(rows,cols,s2) );
00052   m3 = m1;
00053   m3.array() += s2;
00054   VERIFY_IS_APPROX(m3, (m1.array() + s2).matrix());
00055   m3 = m1;
00056   m3.array() -= s1;
00057   VERIFY_IS_APPROX(m3, (m1.array() - s1).matrix());
00058 
00059   // reductions
00060   VERIFY_IS_MUCH_SMALLER_THAN(m1.colwise().sum().sum() - m1.sum(), m1.cwiseAbs().maxCoeff());
00061   VERIFY_IS_MUCH_SMALLER_THAN(m1.rowwise().sum().sum() - m1.sum(), m1.cwiseAbs().maxCoeff());
00062   VERIFY_IS_MUCH_SMALLER_THAN(m1.colwise().sum() + m2.colwise().sum() - (m1+m2).colwise().sum(), (m1+m2).cwiseAbs().maxCoeff());
00063   VERIFY_IS_MUCH_SMALLER_THAN(m1.rowwise().sum() - m2.rowwise().sum() - (m1-m2).rowwise().sum(), (m1-m2).cwiseAbs().maxCoeff());
00064   VERIFY_IS_APPROX(m1.colwise().sum(), m1.colwise().redux(internal::scalar_sum_op<Scalar>()));
00065 
00066   // vector-wise ops
00067   m3 = m1;
00068   VERIFY_IS_APPROX(m3.colwise() += cv1, m1.colwise() + cv1);
00069   m3 = m1;
00070   VERIFY_IS_APPROX(m3.colwise() -= cv1, m1.colwise() - cv1);
00071   m3 = m1;
00072   VERIFY_IS_APPROX(m3.rowwise() += rv1, m1.rowwise() + rv1);
00073   m3 = m1;
00074   VERIFY_IS_APPROX(m3.rowwise() -= rv1, m1.rowwise() - rv1);
00075   
00076   // empty objects
00077   VERIFY_IS_APPROX(m1.block(0,0,0,cols).colwise().sum(),  RowVectorType::Zero(cols));
00078   VERIFY_IS_APPROX(m1.block(0,0,rows,0).rowwise().prod(), ColVectorType::Ones(rows));
00079   
00080   // verify the const accessors exist
00081   const Scalar& ref_m1 = m.matrix().array().coeffRef(0);
00082   const Scalar& ref_m2 = m.matrix().array().coeffRef(0,0);
00083   const Scalar& ref_a1 = m.array().matrix().coeffRef(0);
00084   const Scalar& ref_a2 = m.array().matrix().coeffRef(0,0);
00085   VERIFY(&ref_a1 == &ref_m1);
00086   VERIFY(&ref_a2 == &ref_m2);
00087 }
00088 
00089 template<typename MatrixType> void comparisons(const MatrixType& m)
00090 {
00091   typedef typename MatrixType::Index Index;
00092   typedef typename MatrixType::Scalar Scalar;
00093   typedef typename NumTraits<Scalar>::Real RealScalar;
00094   typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> VectorType;
00095 
00096   Index rows = m.rows();
00097   Index cols = m.cols();
00098 
00099   Index r = internal::random<Index>(0, rows-1),
00100         c = internal::random<Index>(0, cols-1);
00101 
00102   MatrixType m1 = MatrixType::Random(rows, cols),
00103              m2 = MatrixType::Random(rows, cols),
00104              m3(rows, cols);
00105 
00106   VERIFY(((m1.array() + Scalar(1)) > m1.array()).all());
00107   VERIFY(((m1.array() - Scalar(1)) < m1.array()).all());
00108   if (rows*cols>1)
00109   {
00110     m3 = m1;
00111     m3(r,c) += 1;
00112     VERIFY(! (m1.array() < m3.array()).all() );
00113     VERIFY(! (m1.array() > m3.array()).all() );
00114   }
00115 
00116   // comparisons to scalar
00117   VERIFY( (m1.array() != (m1(r,c)+1) ).any() );
00118   VERIFY( (m1.array() > (m1(r,c)-1) ).any() );
00119   VERIFY( (m1.array() < (m1(r,c)+1) ).any() );
00120   VERIFY( (m1.array() == m1(r,c) ).any() );
00121 
00122   // test Select
00123   VERIFY_IS_APPROX( (m1.array()<m2.array()).select(m1,m2), m1.cwiseMin(m2) );
00124   VERIFY_IS_APPROX( (m1.array()>m2.array()).select(m1,m2), m1.cwiseMax(m2) );
00125   Scalar mid = (m1.cwiseAbs().minCoeff() + m1.cwiseAbs().maxCoeff())/Scalar(2);
00126   for (int j=0; j<cols; ++j)
00127   for (int i=0; i<rows; ++i)
00128     m3(i,j) = internal::abs(m1(i,j))<mid ? 0 : m1(i,j);
00129   VERIFY_IS_APPROX( (m1.array().abs()<MatrixType::Constant(rows,cols,mid).array())
00130                         .select(MatrixType::Zero(rows,cols),m1), m3);
00131   // shorter versions:
00132   VERIFY_IS_APPROX( (m1.array().abs()<MatrixType::Constant(rows,cols,mid).array())
00133                         .select(0,m1), m3);
00134   VERIFY_IS_APPROX( (m1.array().abs()>=MatrixType::Constant(rows,cols,mid).array())
00135                         .select(m1,0), m3);
00136   // even shorter version:
00137   VERIFY_IS_APPROX( (m1.array().abs()<mid).select(0,m1), m3);
00138 
00139   // count
00140   VERIFY(((m1.array().abs()+1)>RealScalar(0.1)).count() == rows*cols);
00141 
00142   typedef Matrix<typename MatrixType::Index, Dynamic, 1> VectorOfIndices;
00143 
00144   // TODO allows colwise/rowwise for array
00145   VERIFY_IS_APPROX(((m1.array().abs()+1)>RealScalar(0.1)).matrix().colwise().count(), VectorOfIndices::Constant(cols,rows).transpose());
00146   VERIFY_IS_APPROX(((m1.array().abs()+1)>RealScalar(0.1)).matrix().rowwise().count(), VectorOfIndices::Constant(rows, cols));
00147 }
00148 
00149 template<typename VectorType> void lpNorm(const VectorType& v)
00150 {
00151   VectorType u = VectorType::Random(v.size());
00152 
00153   VERIFY_IS_APPROX(u.template lpNorm<Infinity>(), u.cwiseAbs().maxCoeff());
00154   VERIFY_IS_APPROX(u.template lpNorm<1>(), u.cwiseAbs().sum());
00155   VERIFY_IS_APPROX(u.template lpNorm<2>(), internal::sqrt(u.array().abs().square().sum()));
00156   VERIFY_IS_APPROX(internal::pow(u.template lpNorm<5>(), typename VectorType::RealScalar(5)), u.array().abs().pow(5).sum());
00157 }
00158 
00159 void test_array_for_matrix()
00160 {
00161   int maxsize = 40;
00162   for(int i = 0; i < g_repeat; i++) {
00163     CALL_SUBTEST_1( array_for_matrix(Matrix<float, 1, 1>()) );
00164     CALL_SUBTEST_2( array_for_matrix(Matrix2f()) );
00165     CALL_SUBTEST_3( array_for_matrix(Matrix4d()) );
00166     CALL_SUBTEST_4( array_for_matrix(MatrixXcf(internal::random<int>(1,maxsize), internal::random<int>(1,maxsize))) );
00167     CALL_SUBTEST_5( array_for_matrix(MatrixXf(internal::random<int>(1,maxsize), internal::random<int>(1,maxsize))) );
00168     CALL_SUBTEST_6( array_for_matrix(MatrixXi(internal::random<int>(1,maxsize), internal::random<int>(1,maxsize))) );
00169   }
00170   for(int i = 0; i < g_repeat; i++) {
00171     CALL_SUBTEST_1( comparisons(Matrix<float, 1, 1>()) );
00172     CALL_SUBTEST_2( comparisons(Matrix2f()) );
00173     CALL_SUBTEST_3( comparisons(Matrix4d()) );
00174     CALL_SUBTEST_5( comparisons(MatrixXf(internal::random<int>(1,maxsize), internal::random<int>(1,maxsize))) );
00175     CALL_SUBTEST_6( comparisons(MatrixXi(internal::random<int>(1,maxsize), internal::random<int>(1,maxsize))) );
00176   }
00177   for(int i = 0; i < g_repeat; i++) {
00178     CALL_SUBTEST_1( lpNorm(Matrix<float, 1, 1>()) );
00179     CALL_SUBTEST_2( lpNorm(Vector2f()) );
00180     CALL_SUBTEST_7( lpNorm(Vector3d()) );
00181     CALL_SUBTEST_8( lpNorm(Vector4f()) );
00182     CALL_SUBTEST_5( lpNorm(VectorXf(internal::random<int>(1,maxsize))) );
00183     CALL_SUBTEST_4( lpNorm(VectorXcf(internal::random<int>(1,maxsize))) );
00184   }
00185 }


libicr
Author(s): Robert Krug
autogenerated on Mon Jan 6 2014 11:32:28