prec_inverse_4x4.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) 2009 Benoit Jacob <jacob.benoit.1@gmail.com>
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 #include <Eigen/LU>
00027 #include <algorithm>
00028 
00029 template<typename MatrixType> void inverse_permutation_4x4()
00030 {
00031   typedef typename MatrixType::Scalar Scalar;
00032   typedef typename MatrixType::RealScalar RealScalar;
00033   Vector4i indices(0,1,2,3);
00034   for(int i = 0; i < 24; ++i)
00035   {
00036     MatrixType m = PermutationMatrix<4>(indices);
00037     MatrixType inv = m.inverse();
00038     double error = double( (m*inv-MatrixType::Identity()).norm() / NumTraits<Scalar>::epsilon() );
00039     EIGEN_DEBUG_VAR(error)
00040     VERIFY(error == 0.0);
00041     std::next_permutation(indices.data(),indices.data()+4);
00042   }
00043 }
00044 
00045 template<typename MatrixType> void inverse_general_4x4(int repeat)
00046 {
00047   typedef typename MatrixType::Scalar Scalar;
00048   typedef typename MatrixType::RealScalar RealScalar;
00049   double error_sum = 0., error_max = 0.;
00050   for(int i = 0; i < repeat; ++i)
00051   {
00052     MatrixType m;
00053     RealScalar absdet;
00054     do {
00055       m = MatrixType::Random();
00056       absdet = internal::abs(m.determinant());
00057     } while(absdet < NumTraits<Scalar>::epsilon());
00058     MatrixType inv = m.inverse();
00059     double error = double( (m*inv-MatrixType::Identity()).norm() * absdet / NumTraits<Scalar>::epsilon() );
00060     error_sum += error;
00061     error_max = std::max(error_max, error);
00062   }
00063   std::cerr << "inverse_general_4x4, Scalar = " << type_name<Scalar>() << std::endl;
00064   double error_avg = error_sum / repeat;
00065   EIGEN_DEBUG_VAR(error_avg);
00066   EIGEN_DEBUG_VAR(error_max);
00067    // FIXME that 1.25 used to be a 1.0 until the NumTraits changes on 28 April 2010, what's going wrong??
00068    // FIXME that 1.25 used to be 1.2 until we tested gcc 4.1 on 30 June 2010 and got 1.21.
00069   VERIFY(error_avg < (NumTraits<Scalar>::IsComplex ? 8.0 : 1.25));
00070   VERIFY(error_max < (NumTraits<Scalar>::IsComplex ? 64.0 : 20.0));
00071 }
00072 
00073 void test_prec_inverse_4x4()
00074 {
00075   CALL_SUBTEST_1((inverse_permutation_4x4<Matrix4f>()));
00076   CALL_SUBTEST_1(( inverse_general_4x4<Matrix4f>(200000 * g_repeat) ));
00077 
00078   CALL_SUBTEST_2((inverse_permutation_4x4<Matrix<double,4,4,RowMajor> >()));
00079   CALL_SUBTEST_2(( inverse_general_4x4<Matrix<double,4,4,RowMajor> >(200000 * g_repeat) ));
00080 
00081   CALL_SUBTEST_3((inverse_permutation_4x4<Matrix4cf>()));
00082   CALL_SUBTEST_3((inverse_general_4x4<Matrix4cf>(50000 * g_repeat)));
00083 }


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