aboutsummaryrefslogtreecommitdiffhomepage
path: root/test/determinant.cpp
diff options
context:
space:
mode:
authorGravatar Benoit Jacob <jacob.benoit.1@gmail.com>2009-01-04 15:26:32 +0000
committerGravatar Benoit Jacob <jacob.benoit.1@gmail.com>2009-01-04 15:26:32 +0000
commit15ca6659acea545178116096bff3e42068b4f4cb (patch)
tree0f4093c77a3312b69dff071db7fd1ad14a3d6d25 /test/determinant.cpp
parentd9e5fd393a48db368dd90cf7119ebb3d774111cb (diff)
* the 4th template param of Matrix is now Options. One bit for storage
order, one bit for enabling/disabling auto-alignment. If you want to disable, do: Matrix<float,4,1,Matrix_DontAlign> The Matrix_ prefix is the only way I can see to avoid ambiguity/pollution. The old RowMajor, ColMajor constants are deprecated, remain for now. * this prompted several improvements in matrix_storage. ei_aligned_array renamed to ei_matrix_array and moved there. The %16==0 tests are now much more centralized in 1 place there. * unalignedassert test: updated * update FindEigen2.cmake from KDElibs * determinant test: use VERIFY_IS_APPROX to fix false positives; add testing of 1 big matrix
Diffstat (limited to 'test/determinant.cpp')
-rw-r--r--test/determinant.cpp17
1 files changed, 9 insertions, 8 deletions
diff --git a/test/determinant.cpp b/test/determinant.cpp
index 68e91e4c1..bc647d25d 100644
--- a/test/determinant.cpp
+++ b/test/determinant.cpp
@@ -38,8 +38,8 @@ template<typename MatrixType> void determinant(const MatrixType& m)
m2.setRandom();
typedef typename MatrixType::Scalar Scalar;
Scalar x = ei_random<Scalar>();
- VERIFY(ei_isApprox(MatrixType::Identity(size, size).determinant(), Scalar(1)));
- VERIFY(ei_isApprox((m1*m2).determinant(), m1.determinant() * m2.determinant()));
+ VERIFY_IS_APPROX(MatrixType::Identity(size, size).determinant(), Scalar(1));
+ VERIFY_IS_APPROX((m1*m2).determinant(), m1.determinant() * m2.determinant());
if(size==1) return;
int i = ei_random<int>(0, size-1);
int j;
@@ -48,18 +48,18 @@ template<typename MatrixType> void determinant(const MatrixType& m)
} while(j==i);
m2 = m1;
m2.row(i).swap(m2.row(j));
- VERIFY(ei_isApprox(m2.determinant(), -m1.determinant()));
+ VERIFY_IS_APPROX(m2.determinant(), -m1.determinant());
m2 = m1;
m2.col(i).swap(m2.col(j));
- VERIFY(ei_isApprox(m2.determinant(), -m1.determinant()));
- VERIFY(ei_isApprox(m2.determinant(), m2.transpose().determinant()));
- VERIFY(ei_isApprox(ei_conj(m2.determinant()), m2.adjoint().determinant()));
+ VERIFY_IS_APPROX(m2.determinant(), -m1.determinant());
+ VERIFY_IS_APPROX(m2.determinant(), m2.transpose().determinant());
+ VERIFY_IS_APPROX(ei_conj(m2.determinant()), m2.adjoint().determinant());
m2 = m1;
m2.row(i) += x*m2.row(j);
- VERIFY(ei_isApprox(m2.determinant(), m1.determinant()));
+ VERIFY_IS_APPROX(m2.determinant(), m1.determinant());
m2 = m1;
m2.row(i) *= x;
- VERIFY(ei_isApprox(m2.determinant(), m1.determinant() * x));
+ VERIFY_IS_APPROX(m2.determinant(), m1.determinant() * x);
}
void test_determinant()
@@ -72,4 +72,5 @@ void test_determinant()
CALL_SUBTEST( determinant(Matrix<std::complex<double>, 10, 10>()) );
CALL_SUBTEST( determinant(MatrixXd(20, 20)) );
}
+ CALL_SUBTEST( determinant(MatrixXd(200, 200)) );
}