adjoint.cpp File Reference

#include "main.h"

Classes
struct	adjoint_specific< true >
struct	adjoint_specific< false >

Functions
template<typename MatrixType , typename Scalar = typename MatrixType::Scalar>
MatrixType	RandomMatrix (Index rows, Index cols, Scalar min, Scalar max)
template<typename MatrixType >
void	adjoint (const MatrixType &m)
template<int >
void	adjoint_extra ()
	EIGEN_DECLARE_TEST (adjoint)

Function Documentation

adjoint()

template<typename MatrixType >

void adjoint

(

const MatrixType &

m

)

                                  {
  /* this test covers the following files:
     Transpose.h Conjugate.h Dot.h
  */
  using std::abs;
  typedef typename MatrixType::Scalar Scalar;
  typedef typename NumTraits<Scalar>::Real RealScalar;
  typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> VectorType;
  typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, MatrixType::RowsAtCompileTime> SquareMatrixType;
  const Index PacketSize = internal::packet_traits<Scalar>::size;
 
  Index rows = m.rows();
  Index cols = m.cols();
 
  // Avoid integer overflow by limiting input values.
  RealScalar rmin = static_cast<RealScalar>(NumTraits<Scalar>::IsInteger ? NumTraits<Scalar>::IsSigned ? -100 : 0 : -1);
  RealScalar rmax = static_cast<RealScalar>(NumTraits<Scalar>::IsInteger ? 100 : 1);
 
  MatrixType m1 = RandomMatrix<MatrixType>(rows, cols, rmin, rmax),
             m2 = RandomMatrix<MatrixType>(rows, cols, rmin, rmax), m3(rows, cols),
             square = RandomMatrix<SquareMatrixType>(rows, rows, rmin, rmax);
  VectorType v1 = RandomMatrix<VectorType>(rows, 1, rmin, rmax), v2 = RandomMatrix<VectorType>(rows, 1, rmin, rmax),
             v3 = RandomMatrix<VectorType>(rows, 1, rmin, rmax), vzero = VectorType::Zero(rows);
 
  Scalar s1 = internal::random<Scalar>(rmin, rmax), s2 = internal::random<Scalar>(rmin, rmax);
 
  // check basic compatibility of adjoint, transpose, conjugate
  VERIFY_IS_APPROX(m1.transpose().conjugate().adjoint(), m1);
  VERIFY_IS_APPROX(m1.adjoint().conjugate().transpose(), m1);
 
  // check multiplicative behavior
  VERIFY_IS_APPROX((m1.adjoint() * m2).adjoint(), m2.adjoint() * m1);
  VERIFY_IS_APPROX((s1 * m1).adjoint(), numext::conj(s1) * m1.adjoint());
 
  // check basic properties of dot, squaredNorm
  VERIFY_IS_APPROX(numext::conj(v1.dot(v2)), v2.dot(v1));
  VERIFY_IS_APPROX(numext::real(v1.dot(v1)), v1.squaredNorm());
 
  adjoint_specific<NumTraits<Scalar>::IsInteger>::run(v1, v2, v3, square, s1, s2);
 
  VERIFY_IS_MUCH_SMALLER_THAN(abs(vzero.dot(v1)), static_cast<RealScalar>(1));
 
  // like in testBasicStuff, test operator() to check const-qualification
  Index r = internal::random<Index>(0, rows - 1), c = internal::random<Index>(0, cols - 1);
  VERIFY_IS_APPROX(m1.conjugate()(r, c), numext::conj(m1(r, c)));
  VERIFY_IS_APPROX(m1.adjoint()(c, r), numext::conj(m1(r, c)));
 
  // check inplace transpose
  m3 = m1;
  m3.transposeInPlace();
  VERIFY_IS_APPROX(m3, m1.transpose());
  m3.transposeInPlace();
  VERIFY_IS_APPROX(m3, m1);
 
  if (PacketSize < m3.rows() && PacketSize < m3.cols()) {
    m3 = m1;
    Index i = internal::random<Index>(0, m3.rows() - PacketSize);
    Index j = internal::random<Index>(0, m3.cols() - PacketSize);
    m3.template block<PacketSize, PacketSize>(i, j).transposeInPlace();
    VERIFY_IS_APPROX((m3.template block<PacketSize, PacketSize>(i, j)),
                     (m1.template block<PacketSize, PacketSize>(i, j).transpose()));
    m3.template block<PacketSize, PacketSize>(i, j).transposeInPlace();
    VERIFY_IS_APPROX(m3, m1);
  }
 
  // check inplace adjoint
  m3 = m1;
  m3.adjointInPlace();
  VERIFY_IS_APPROX(m3, m1.adjoint());
  m3.transposeInPlace();
  VERIFY_IS_APPROX(m3, m1.conjugate());
 
  // check mixed dot product
  typedef Matrix<RealScalar, MatrixType::RowsAtCompileTime, 1> RealVectorType;
  RealVectorType rv1 = RandomMatrix<RealVectorType>(rows, 1, rmin, rmax);
 
  VERIFY_IS_APPROX(v1.dot(rv1.template cast<Scalar>()), v1.dot(rv1));
  VERIFY_IS_APPROX(rv1.template cast<Scalar>().dot(v1), rv1.dot(v1));
 
  VERIFY(is_same_type(m1, m1.template conjugateIf<false>()));
  VERIFY(is_same_type(m1.conjugate(), m1.template conjugateIf<true>()));
}

References abs(), calibrate::c, cols, conj(), i, Eigen::is_same_type(), j, m, m1, m2(), UniformPSDSelfTest::r, Global_Parameters::rmax, Global_Parameters::rmin, rows, run(), size, Eigen::square(), v1(), v2(), VERIFY, VERIFY_IS_APPROX, VERIFY_IS_MUCH_SMALLER_THAN, and oomph::PseudoSolidHelper::Zero.

Referenced by Eigen::CompleteOrthogonalDecomposition< MatrixType_, PermutationIndex_ >::_solve_impl(), Eigen::MatrixBase< Derived >::adjointInPlace(), Eigen::CompleteOrthogonalDecomposition< MatrixType_, PermutationIndex_ >::applyZAdjointOnTheLeftInPlace(), block(), Eigen::DGMRES< MatrixType_, Preconditioner_ >::dgmresCycle(), EIGEN_DECLARE_TEST(), Eigen::internal::gmres(), Eigen::internal::idrs(), Eigen::internal::idrstabl(), product_extra(), product_selfadjoint(), Eigen::internal::solve_assertion< CwiseUnaryOp< Eigen::internal::scalar_conjugate_op< Scalar >, const Transpose< Derived > > >::run(), symm(), and trmm().

adjoint_extra()

template<int >

void adjoint_extra

(

)

                     {
  MatrixXcf a(10, 10), b(10, 10);
  VERIFY_RAISES_ASSERT(a = a.transpose());
  VERIFY_RAISES_ASSERT(a = a.transpose() + b);
  VERIFY_RAISES_ASSERT(a = b + a.transpose());
  VERIFY_RAISES_ASSERT(a = a.conjugate().transpose());
  VERIFY_RAISES_ASSERT(a = a.adjoint());
  VERIFY_RAISES_ASSERT(a = a.adjoint() + b);
  VERIFY_RAISES_ASSERT(a = b + a.adjoint());
 
  // no assertion should be triggered for these cases:
  a.transpose() = a.transpose();
  a.transpose() += a.transpose();
  a.transpose() += a.transpose() + b;
  a.transpose() = a.adjoint();
  a.transpose() += a.adjoint();
  a.transpose() += a.adjoint() + b;
 
  // regression tests for check_for_aliasing
  MatrixXd c(10, 10);
  c = 1.0 * MatrixXd::Ones(10, 10) + c;
  c = MatrixXd::Ones(10, 10) * 1.0 + c;
  c = c + MatrixXd::Ones(10, 10).cwiseProduct(MatrixXd::Zero(10, 10));
  c = MatrixXd::Ones(10, 10) * MatrixXd::Zero(10, 10);
 
  // regression for bug 1646
  for (int j = 0; j < 10; ++j) {
    c.col(j).head(j) = c.row(j).head(j);
  }
 
  for (int j = 0; j < 10; ++j) {
    c.col(j) = c.row(j);
  }
 
  a.conservativeResize(1, 1);
  a = a.transpose();
 
  a.conservativeResize(0, 0);
  a = a.transpose();
}

References a, b, calibrate::c, j, VERIFY_RAISES_ASSERT, and oomph::PseudoSolidHelper::Zero.

EIGEN_DECLARE_TEST()

EIGEN_DECLARE_TEST

(

adjoint

)

                            {
  for (int i = 0; i < g_repeat; i++) {
    CALL_SUBTEST_1(adjoint(Matrix<float, 1, 1>()));
    CALL_SUBTEST_2(adjoint(Matrix3d()));
    CALL_SUBTEST_3(adjoint(Matrix4f()));
 
    CALL_SUBTEST_4(adjoint(MatrixXcf(internal::random<int>(1, EIGEN_TEST_MAX_SIZE / 2),
                                     internal::random<int>(1, EIGEN_TEST_MAX_SIZE / 2))));
    CALL_SUBTEST_5(adjoint(
        MatrixXi(internal::random<int>(1, EIGEN_TEST_MAX_SIZE), internal::random<int>(1, EIGEN_TEST_MAX_SIZE))));
    CALL_SUBTEST_6(adjoint(
        MatrixXf(internal::random<int>(1, EIGEN_TEST_MAX_SIZE), internal::random<int>(1, EIGEN_TEST_MAX_SIZE))));
 
    // Complement for 128 bits vectorization:
    CALL_SUBTEST_8(adjoint(Matrix2d()));
    CALL_SUBTEST_9(adjoint(Matrix<int, 4, 4>()));
 
    // 256 bits vectorization:
    CALL_SUBTEST_10(adjoint(Matrix<float, 8, 8>()));
    CALL_SUBTEST_11(adjoint(Matrix<double, 4, 4>()));
    CALL_SUBTEST_12(adjoint(Matrix<int, 8, 8>()));
  }
  // test a large static matrix only once
  CALL_SUBTEST_7(adjoint(Matrix<float, 100, 100>()));
 
  CALL_SUBTEST_13(adjoint_extra<0>());
}

References adjoint(), CALL_SUBTEST_1, CALL_SUBTEST_10, CALL_SUBTEST_11, CALL_SUBTEST_12, CALL_SUBTEST_13, CALL_SUBTEST_2, CALL_SUBTEST_3, CALL_SUBTEST_4, CALL_SUBTEST_5, CALL_SUBTEST_6, CALL_SUBTEST_7, CALL_SUBTEST_8, CALL_SUBTEST_9, EIGEN_TEST_MAX_SIZE, Eigen::g_repeat, and i.

RandomMatrix()

template<typename MatrixType , typename Scalar = typename MatrixType::Scalar>

MatrixType RandomMatrix	(	Index	rows,
		Index	cols,
		Scalar	min,
		Scalar	max
	)

                                                                        {
  MatrixType M = MatrixType(rows, cols);
  for (Index i = 0; i < rows; ++i) {
    for (Index j = 0; j < cols; ++j) {
      M(i, j) = Eigen::internal::random<Scalar>(min, max);
    }
  }
  return M;
}

References cols, i, j, max, min, and rows.