C++ VERIFY_IS_APPROX函数代码示例

template<typename VectorType> void vectorVisitor(const VectorType& w)
{
  typedef typename VectorType::Scalar Scalar;
  typedef typename VectorType::Index Index;

  Index size = w.size();

  // construct a random vector where all coefficients are different
  VectorType v;
  v = VectorType::Random(size);
  for(Index i = 0; i < size; i++)
    for(Index i2 = 0; i2 < i; i2++)
      while(v(i) == v(i2)) // yes, ==
        v(i) = ei_random<Scalar>();
  
  Scalar minc = Scalar(1000), maxc = Scalar(-1000);
  Index minidx=0,maxidx=0;
  for(Index i = 0; i < size; i++)
  {
    if(v(i) < minc)
    {
      minc = v(i);
      minidx = i;
    }
    if(v(i) > maxc)
    {
      maxc = v(i);
      maxidx = i;
    }
  }
  Index eigen_minidx, eigen_maxidx;
  Scalar eigen_minc, eigen_maxc;
  eigen_minc = v.minCoeff(&eigen_minidx);
  eigen_maxc = v.maxCoeff(&eigen_maxidx);
  VERIFY(minidx == eigen_minidx);
  VERIFY(maxidx == eigen_maxidx);
  VERIFY_IS_APPROX(minc, eigen_minc);
  VERIFY_IS_APPROX(maxc, eigen_maxc);
  VERIFY_IS_APPROX(minc, v.minCoeff());
  VERIFY_IS_APPROX(maxc, v.maxCoeff());
}

template<typename Scalar> void smallVectors()
{
  typedef Matrix<Scalar, 1, 2> V2;
  typedef Matrix<Scalar, 3, 1> V3;
  typedef Matrix<Scalar, 1, 4> V4;
  typedef Matrix<Scalar, Dynamic, 1> VX;
  Scalar x1 = internal::random<Scalar>(),
         x2 = internal::random<Scalar>(),
         x3 = internal::random<Scalar>(),
         x4 = internal::random<Scalar>();
  V2 v2(x1, x2);
  V3 v3(x1, x2, x3);
  V4 v4(x1, x2, x3, x4);
  VERIFY_IS_APPROX(x1, v2.x());
  VERIFY_IS_APPROX(x1, v3.x());
  VERIFY_IS_APPROX(x1, v4.x());
  VERIFY_IS_APPROX(x2, v2.y());
  VERIFY_IS_APPROX(x2, v3.y());
  VERIFY_IS_APPROX(x2, v4.y());
  VERIFY_IS_APPROX(x3, v3.z());
  VERIFY_IS_APPROX(x3, v4.z());
  VERIFY_IS_APPROX(x4, v4.w());

  if (!NumTraits<Scalar>::IsInteger)
  {
    VERIFY_RAISES_ASSERT(V3(2, 1))
    VERIFY_RAISES_ASSERT(V3(3, 2))
    VERIFY_RAISES_ASSERT(V3(Scalar(3), 1))
    VERIFY_RAISES_ASSERT(V3(3, Scalar(1)))
    VERIFY_RAISES_ASSERT(V3(Scalar(3), Scalar(1)))
    VERIFY_RAISES_ASSERT(V3(Scalar(123), Scalar(123)))

    VERIFY_RAISES_ASSERT(V4(1, 3))
    VERIFY_RAISES_ASSERT(V4(2, 4))
    VERIFY_RAISES_ASSERT(V4(1, Scalar(4)))
    VERIFY_RAISES_ASSERT(V4(Scalar(1), 4))
    VERIFY_RAISES_ASSERT(V4(Scalar(1), Scalar(4)))
    VERIFY_RAISES_ASSERT(V4(Scalar(123), Scalar(123)))

    VERIFY_RAISES_ASSERT(VX(3, 2))
    VERIFY_RAISES_ASSERT(VX(Scalar(3), 1))
    VERIFY_RAISES_ASSERT(VX(3, Scalar(1)))
    VERIFY_RAISES_ASSERT(VX(Scalar(3), Scalar(1)))
    VERIFY_RAISES_ASSERT(VX(Scalar(123), Scalar(123)))
  }
}

template<typename MatrixType> void swap(const MatrixType& m)
{
  typedef typename other_matrix_type<MatrixType>::type OtherMatrixType;
  typedef typename MatrixType::Scalar Scalar;

  eigen_assert((!internal::is_same<MatrixType,OtherMatrixType>::value));
  typename MatrixType::Index rows = m.rows();
  typename MatrixType::Index cols = m.cols();

  // construct 3 matrix guaranteed to be distinct
  MatrixType m1 = MatrixType::Random(rows,cols);
  MatrixType m2 = MatrixType::Random(rows,cols) + Scalar(100) * MatrixType::Identity(rows,cols);
  OtherMatrixType m3 = OtherMatrixType::Random(rows,cols) + Scalar(200) * OtherMatrixType::Identity(rows,cols);

  MatrixType m1_copy = m1;
  MatrixType m2_copy = m2;
  OtherMatrixType m3_copy = m3;

  // test swapping 2 matrices of same type
  m1.swap(m2);
  VERIFY_IS_APPROX(m1,m2_copy);
  VERIFY_IS_APPROX(m2,m1_copy);
  m1 = m1_copy;
  m2 = m2_copy;

  // test swapping 2 matrices of different types
  m1.swap(m3);
  VERIFY_IS_APPROX(m1,m3_copy);
  VERIFY_IS_APPROX(m3,m1_copy);
  m1 = m1_copy;
  m3 = m3_copy;

  // test swapping matrix with expression
  m1.swap(m2.block(0,0,rows,cols));
  VERIFY_IS_APPROX(m1,m2_copy);
  VERIFY_IS_APPROX(m2,m1_copy);
  m1 = m1_copy;
  m2 = m2_copy;

  // test swapping two expressions of different types
  m1.transpose().swap(m3.transpose());
  VERIFY_IS_APPROX(m1,m3_copy);
  VERIFY_IS_APPROX(m3,m1_copy);
  m1 = m1_copy;
  m3 = m3_copy;

  // test assertion on mismatching size -- matrix case
  VERIFY_RAISES_ASSERT(m1.swap(m1.row(0)));
  // test assertion on mismatching size -- xpr case
  VERIFY_RAISES_ASSERT(m1.row(0).swap(m1));
}

void verify_is_approx_upto_permutation(const VectorType& vec1, const VectorType& vec2)
{
  typedef typename NumTraits<typename VectorType::Scalar>::Real RealScalar;

  VERIFY(vec1.cols() == 1);
  VERIFY(vec2.cols() == 1);
  VERIFY(vec1.rows() == vec2.rows());
  for (int k = 1; k <= vec1.rows(); ++k)
  {
    VERIFY_IS_APPROX(vec1.array().pow(RealScalar(k)).sum(), vec2.array().pow(RealScalar(k)).sum());
  }
}

template<typename MatrixType> void selfadjoint(const MatrixType& m)
{
  typedef typename MatrixType::Index Index;
  typedef typename MatrixType::Scalar Scalar;

  Index rows = m.rows();
  Index cols = m.cols();

  MatrixType m1 = MatrixType::Random(rows, cols),
             m2 = MatrixType::Random(rows, cols),
             m3(rows, cols),
             m4(rows, cols);

  m1.diagonal() = m1.diagonal().real().template cast<Scalar>();

  // check selfadjoint to dense
  m3 = m1.template selfadjointView<Upper>();
  VERIFY_IS_APPROX(MatrixType(m3.template triangularView<Upper>()), MatrixType(m1.template triangularView<Upper>()));
  VERIFY_IS_APPROX(m3, m3.adjoint());

  m3 = m1.template selfadjointView<Lower>();
  VERIFY_IS_APPROX(MatrixType(m3.template triangularView<Lower>()), MatrixType(m1.template triangularView<Lower>()));
  VERIFY_IS_APPROX(m3, m3.adjoint());

  m3 = m1.template selfadjointView<Upper>();
  m4 = m2;
  m4 += m1.template selfadjointView<Upper>();
  VERIFY_IS_APPROX(m4, m2+m3);

  m3 = m1.template selfadjointView<Lower>();
  m4 = m2;
  m4 -= m1.template selfadjointView<Lower>();
  VERIFY_IS_APPROX(m4, m2-m3);
}

template<typename Func> void adolc_forward_jacobian(const Func& f)
{
    typename Func::InputType x = Func::InputType::Random(f.inputs());
    typename Func::ValueType y(f.values()), yref(f.values());
    typename Func::JacobianType j(f.values(),f.inputs()), jref(f.values(),f.inputs());

    jref.setZero();
    yref.setZero();
    f(x,&yref,&jref);
//     std::cerr << y.transpose() << "\n\n";;
//     std::cerr << j << "\n\n";;

    j.setZero();
    y.setZero();
    AdolcForwardJacobian<Func> autoj(f);
    autoj(x, &y, &j);
//     std::cerr << y.transpose() << "\n\n";;
//     std::cerr << j << "\n\n";;

    VERIFY_IS_APPROX(y, yref);
    VERIFY_IS_APPROX(j, jref);
}

void check_stdlist_matrix(const MatrixType& m)
{
  typename MatrixType::Index rows = m.rows();
  typename MatrixType::Index cols = m.cols();
  MatrixType x = MatrixType::Random(rows,cols), y = MatrixType::Random(rows,cols);
  std::list<MatrixType> v(10, MatrixType(rows,cols)), w(20, y);
  typename std::list<MatrixType>::iterator itv = get(v, 5);
  typename std::list<MatrixType>::iterator itw = get(w, 6);
  *itv = x;
  *itw = *itv;
  VERIFY_IS_APPROX(*itw, *itv);
  v = w;
  itv = v.begin();
  itw = w.begin();
  for(int i = 0; i < 20; i++)
  {
    VERIFY_IS_APPROX(*itw, *itv);
    ++itv;
    ++itw;
  }

  v.resize(21);
  set(v, 20, x);
  VERIFY_IS_APPROX(*get(v, 20), x);
  v.resize(22,y);
  VERIFY_IS_APPROX(*get(v, 21), y);
  v.push_back(x);
  VERIFY_IS_APPROX(*get(v, 22), x);

  // do a lot of push_back such that the list gets internally resized
  // (with memory reallocation)
  MatrixType* ref = &(*get(w, 0));
  for(int i=0; i<30 || ((ref==&(*get(w, 0))) && i<300); ++i)
    v.push_back(*get(w, i%w.size()));
  for(unsigned int i=23; i<v.size(); ++i)
  {
    VERIFY((*get(v, i))==(*get(w, (i-23)%w.size())));
  }
}

template<typename Scalar> void mapQuaternion(void){
  typedef Map<Quaternion<Scalar>, Aligned> MQuaternionA;
  typedef Map<Quaternion<Scalar> > MQuaternionUA;
  typedef Quaternion<Scalar> Quaternionx;

	EIGEN_ALIGN16 Scalar array1[4];
	EIGEN_ALIGN16 Scalar array2[4];
	EIGEN_ALIGN16 Scalar array3[4+1];
	Scalar* array3unaligned = array3+1;

  MQuaternionA(array1).coeffs().setRandom();
  (MQuaternionA(array2)) = MQuaternionA(array1);
  (MQuaternionUA(array3unaligned)) = MQuaternionA(array1);

  Quaternionx q1 = MQuaternionA(array1);
  Quaternionx q2 = MQuaternionA(array2);
  Quaternionx q3 = MQuaternionUA(array3unaligned);

  VERIFY_IS_APPROX(q1.coeffs(), q2.coeffs());
  VERIFY_IS_APPROX(q1.coeffs(), q3.coeffs());
  VERIFY_RAISES_ASSERT((MQuaternionA(array3unaligned)));
}

template <typename MatrixType> void run_nesting_ops(const MatrixType& _m)
{
  typename MatrixType::Nested m(_m);
  typedef typename MatrixType::Scalar Scalar;

#ifdef NDEBUG
  const bool is_debug = false;
#else
  const bool is_debug = true;
#endif

  // Make really sure that we are in debug mode! We don't want any type of
  // inlining for these tests to pass.
  VERIFY(is_debug);

  // The only intention of these tests is to ensure that this code does
  // not trigger any asserts or segmentation faults... more to come.
  VERIFY_IS_APPROX( (m.transpose() * m).diagonal().sum(), (m.transpose() * m).diagonal().sum() );
  VERIFY_IS_APPROX( (m.transpose() * m).diagonal().array().abs().sum(), (m.transpose() * m).diagonal().array().abs().sum() );

  VERIFY_IS_APPROX( (m.transpose() * m).array().abs().sum(), (m.transpose() * m).array().abs().sum() );
}

void check_stdlist_transform(const TransformType&)
{
  typedef typename TransformType::MatrixType MatrixType;
  TransformType x(MatrixType::Random()), y(MatrixType::Random());
  std::list<TransformType> v(10), w(20, y);
  typename std::list<TransformType>::iterator itv = get(v, 5);
  typename std::list<TransformType>::iterator itw = get(w, 6);
  *itv = x;
  *itw = *itv;
  VERIFY_IS_APPROX(*itw, *itv);
  v = w;
  itv = v.begin();
  itw = w.begin();
  for(int i = 0; i < 20; i++)
  {
    VERIFY_IS_APPROX(*itw, *itv);
    ++itv;
    ++itw;
  }

  v.resize(21);
  set(v, 20, x);
  VERIFY_IS_APPROX(*get(v, 20), x);
  v.resize(22,y);
  VERIFY_IS_APPROX(*get(v, 21), y);
  v.push_back(x);
  VERIFY_IS_APPROX(*get(v, 22), x);

  // do a lot of push_back such that the list gets internally resized
  // (with memory reallocation)
  TransformType* ref = &(*get(w, 0));
  for(int i=0; i<30 || ((ref==&(*get(w, 0))) && i<300); ++i)
    v.push_back(*get(w, i%w.size()));
  for(unsigned int i=23; i<v.size(); ++i)
  {
    VERIFY(get(v, i)->matrix()==get(w, (i-23)%w.size())->matrix());
  }
}

template<typename MatrixType> void schur(int size = MatrixType::ColsAtCompileTime)
{
  typedef typename ComplexSchur<MatrixType>::ComplexScalar ComplexScalar;
  typedef typename ComplexSchur<MatrixType>::ComplexMatrixType ComplexMatrixType;

  // Test basic functionality: T is triangular and A = U T U*
  for(int counter = 0; counter < g_repeat; ++counter) {
    MatrixType A = MatrixType::Random(size, size);
    ComplexSchur<MatrixType> schurOfA(A);
    VERIFY_IS_EQUAL(schurOfA.info(), Success);
    ComplexMatrixType U = schurOfA.matrixU();
    ComplexMatrixType T = schurOfA.matrixT();
    for(int row = 1; row < size; ++row) {
      for(int col = 0; col < row; ++col) {
	VERIFY(T(row,col) == (typename MatrixType::Scalar)0);
      }
    }
    VERIFY_IS_APPROX(A.template cast<ComplexScalar>(), U * T * U.adjoint());
  }

  // Test asserts when not initialized
  ComplexSchur<MatrixType> csUninitialized;
  VERIFY_RAISES_ASSERT(csUninitialized.matrixT());
  VERIFY_RAISES_ASSERT(csUninitialized.matrixU());
  VERIFY_RAISES_ASSERT(csUninitialized.info());
  
  // Test whether compute() and constructor returns same result
  MatrixType A = MatrixType::Random(size, size);
  ComplexSchur<MatrixType> cs1;
  cs1.compute(A);
  ComplexSchur<MatrixType> cs2(A);
  VERIFY_IS_EQUAL(cs1.info(), Success);
  VERIFY_IS_EQUAL(cs2.info(), Success);
  VERIFY_IS_EQUAL(cs1.matrixT(), cs2.matrixT());
  VERIFY_IS_EQUAL(cs1.matrixU(), cs2.matrixU());

  // Test computation of only T, not U
  ComplexSchur<MatrixType> csOnlyT(A, false);
  VERIFY_IS_EQUAL(csOnlyT.info(), Success);
  VERIFY_IS_EQUAL(cs1.matrixT(), csOnlyT.matrixT());
  VERIFY_RAISES_ASSERT(csOnlyT.matrixU());

  if (size > 1)
  {
    // Test matrix with NaN
    A(0,0) = std::numeric_limits<typename MatrixType::RealScalar>::quiet_NaN();
    ComplexSchur<MatrixType> csNaN(A);
    VERIFY_IS_EQUAL(csNaN.info(), NoConvergence);
  }
}

template<typename Scalar> void quaternionAlignment(void){
  typedef Quaternion<Scalar,AutoAlign> QuaternionA;
  typedef Quaternion<Scalar,DontAlign> QuaternionUA;

  EIGEN_ALIGN16 Scalar array1[4];
  EIGEN_ALIGN16 Scalar array2[4];
  EIGEN_ALIGN16 Scalar array3[4+1];
  Scalar* arrayunaligned = array3+1;

  QuaternionA *q1 = ::new(reinterpret_cast<void*>(array1)) QuaternionA;
  QuaternionUA *q2 = ::new(reinterpret_cast<void*>(array2)) QuaternionUA;
  QuaternionUA *q3 = ::new(reinterpret_cast<void*>(arrayunaligned)) QuaternionUA;

  q1->coeffs().setRandom();
  *q2 = *q1;
  *q3 = *q1;

  VERIFY_IS_APPROX(q1->coeffs(), q2->coeffs());
  VERIFY_IS_APPROX(q1->coeffs(), q3->coeffs());
  #ifdef EIGEN_VECTORIZE
  VERIFY_RAISES_ASSERT((::new(reinterpret_cast<void*>(arrayunaligned)) QuaternionA));
  #endif
}

template<typename MatrixType> void basicStuffComplex(const MatrixType& m)
{
  typedef typename MatrixType::Index Index;
  typedef typename MatrixType::Scalar Scalar;
  typedef typename NumTraits<Scalar>::Real RealScalar;
  typedef Matrix<RealScalar, MatrixType::RowsAtCompileTime, MatrixType::ColsAtCompileTime> RealMatrixType;

  Index rows = m.rows();
  Index cols = m.cols();

  Scalar s1 = internal::random<Scalar>(),
         s2 = internal::random<Scalar>();

  VERIFY(numext::real(s1)==numext::real_ref(s1));
  VERIFY(numext::imag(s1)==numext::imag_ref(s1));
  numext::real_ref(s1) = numext::real(s2);
  numext::imag_ref(s1) = numext::imag(s2);
  VERIFY(internal::isApprox(s1, s2, NumTraits<RealScalar>::epsilon()));
  // extended precision in Intel FPUs means that s1 == s2 in the line above is not guaranteed.

  RealMatrixType rm1 = RealMatrixType::Random(rows,cols),
                 rm2 = RealMatrixType::Random(rows,cols);
  MatrixType cm(rows,cols);
  cm.real() = rm1;
  cm.imag() = rm2;
  VERIFY_IS_APPROX(static_cast<const MatrixType&>(cm).real(), rm1);
  VERIFY_IS_APPROX(static_cast<const MatrixType&>(cm).imag(), rm2);
  rm1.setZero();
  rm2.setZero();
  rm1 = cm.real();
  rm2 = cm.imag();
  VERIFY_IS_APPROX(static_cast<const MatrixType&>(cm).real(), rm1);
  VERIFY_IS_APPROX(static_cast<const MatrixType&>(cm).imag(), rm2);
  cm.real().setZero();
  VERIFY(static_cast<const MatrixType&>(cm).real().isZero());
  VERIFY(!static_cast<const MatrixType&>(cm).imag().isZero());
}

template<typename BoxType> void alignedbox(const BoxType& _box)
{
  /* this test covers the following files:
     AlignedBox.h
  */

  const int dim = _box.dim();
  typedef typename BoxType::Scalar Scalar;
  typedef typename NumTraits<Scalar>::Real RealScalar;
  typedef Matrix<Scalar, BoxType::AmbientDimAtCompileTime, 1> VectorType;

  VectorType p0 = VectorType::Random(dim);
  VectorType p1 = VectorType::Random(dim);
  RealScalar s1 = ei_random<RealScalar>(0,1);

  BoxType b0(dim);
  BoxType b1(VectorType::Random(dim),VectorType::Random(dim));
  BoxType b2;

  b0.extend(p0);
  b0.extend(p1);
  VERIFY(b0.contains(p0*s1+(Scalar(1)-s1)*p1));
  VERIFY(!b0.contains(p0 + (1+s1)*(p1-p0)));

  (b2 = b0).extend(b1);
  VERIFY(b2.contains(b0));
  VERIFY(b2.contains(b1));
  VERIFY_IS_APPROX(b2.clamp(b0), b0);

  // casting
  const int Dim = BoxType::AmbientDimAtCompileTime;
  typedef typename GetDifferentType<Scalar>::type OtherScalar;
  AlignedBox<OtherScalar,Dim> hp1f = b0.template cast<OtherScalar>();
  VERIFY_IS_APPROX(hp1f.template cast<Scalar>(),b0);
  AlignedBox<Scalar,Dim> hp1d = b0.template cast<Scalar>();
  VERIFY_IS_APPROX(hp1d.template cast<Scalar>(),b0);
}

template<typename LineType> void parametrizedline(const LineType& _line)
{
  /* this test covers the following files:
     ParametrizedLine.h
  */
  typedef typename LineType::Index Index;
  const Index dim = _line.dim();
  typedef typename LineType::Scalar Scalar;
  typedef typename NumTraits<Scalar>::Real RealScalar;
  typedef Matrix<Scalar, LineType::AmbientDimAtCompileTime, 1> VectorType;
  typedef Matrix<Scalar, LineType::AmbientDimAtCompileTime,
                         LineType::AmbientDimAtCompileTime> MatrixType;

  VectorType p0 = VectorType::Random(dim);
  VectorType p1 = VectorType::Random(dim);

  VectorType d0 = VectorType::Random(dim).normalized();

  LineType l0(p0, d0);

  Scalar s0 = internal::random<Scalar>();
  Scalar s1 = internal::abs(internal::random<Scalar>());

  VERIFY_IS_MUCH_SMALLER_THAN( l0.distance(p0), RealScalar(1) );
  VERIFY_IS_MUCH_SMALLER_THAN( l0.distance(p0+s0*d0), RealScalar(1) );
  VERIFY_IS_APPROX( (l0.projection(p1)-p1).norm(), l0.distance(p1) );
  VERIFY_IS_MUCH_SMALLER_THAN( l0.distance(l0.projection(p1)), RealScalar(1) );
  VERIFY_IS_APPROX( Scalar(l0.distance((p0+s0*d0) + d0.unitOrthogonal() * s1)), s1 );

  // casting
  const int Dim = LineType::AmbientDimAtCompileTime;
  typedef typename GetDifferentType<Scalar>::type OtherScalar;
  ParametrizedLine<OtherScalar,Dim> hp1f = l0.template cast<OtherScalar>();
  VERIFY_IS_APPROX(hp1f.template cast<Scalar>(),l0);
  ParametrizedLine<Scalar,Dim> hp1d = l0.template cast<Scalar>();
  VERIFY_IS_APPROX(hp1d.template cast<Scalar>(),l0);
}