1// This file is part of Eigen, a lightweight C++ template library 2// for linear algebra. 3// 4// Copyright (C) 2016 Gael Guennebaud <gael.guennebaud@inria.fr> 5// 6// This Source Code Form is subject to the terms of the Mozilla 7// Public License v. 2.0. If a copy of the MPL was not distributed 8// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 9 10#include <sstream> 11 12#ifdef EIGEN_TEST_MAX_SIZE 13#undef EIGEN_TEST_MAX_SIZE 14#endif 15 16#define EIGEN_TEST_MAX_SIZE 50 17 18#ifdef EIGEN_TEST_PART_1 19#include "cholesky.cpp" 20#endif 21 22#ifdef EIGEN_TEST_PART_2 23#include "lu.cpp" 24#endif 25 26#ifdef EIGEN_TEST_PART_3 27#include "qr.cpp" 28#endif 29 30#ifdef EIGEN_TEST_PART_4 31#include "qr_colpivoting.cpp" 32#endif 33 34#ifdef EIGEN_TEST_PART_5 35#include "qr_fullpivoting.cpp" 36#endif 37 38#ifdef EIGEN_TEST_PART_6 39#include "eigensolver_selfadjoint.cpp" 40#endif 41 42#ifdef EIGEN_TEST_PART_7 43#include "eigensolver_generic.cpp" 44#endif 45 46#ifdef EIGEN_TEST_PART_8 47#include "eigensolver_generalized_real.cpp" 48#endif 49 50#ifdef EIGEN_TEST_PART_9 51#include "jacobisvd.cpp" 52#endif 53 54#ifdef EIGEN_TEST_PART_10 55#include "bdcsvd.cpp" 56#endif 57 58#include <Eigen/Dense> 59 60#undef min 61#undef max 62#undef isnan 63#undef isinf 64#undef isfinite 65 66#include <boost/multiprecision/cpp_dec_float.hpp> 67#include <boost/multiprecision/number.hpp> 68#include <boost/math/special_functions.hpp> 69#include <boost/math/complex.hpp> 70 71namespace mp = boost::multiprecision; 72typedef mp::number<mp::cpp_dec_float<100>, mp::et_on> Real; 73 74namespace Eigen { 75 template<> struct NumTraits<Real> : GenericNumTraits<Real> { 76 static inline Real dummy_precision() { return 1e-50; } 77 }; 78 79 template<typename T1,typename T2,typename T3,typename T4,typename T5> 80 struct NumTraits<boost::multiprecision::detail::expression<T1,T2,T3,T4,T5> > : NumTraits<Real> {}; 81 82 template<> 83 Real test_precision<Real>() { return 1e-50; } 84 85 // needed in C++93 mode where number does not support explicit cast. 86 namespace internal { 87 template<typename NewType> 88 struct cast_impl<Real,NewType> { 89 static inline NewType run(const Real& x) { 90 return x.template convert_to<NewType>(); 91 } 92 }; 93 94 template<> 95 struct cast_impl<Real,std::complex<Real> > { 96 static inline std::complex<Real> run(const Real& x) { 97 return std::complex<Real>(x); 98 } 99 }; 100 } 101} 102 103namespace boost { 104namespace multiprecision { 105 // to make ADL works as expected: 106 using boost::math::isfinite; 107 using boost::math::isnan; 108 using boost::math::isinf; 109 using boost::math::copysign; 110 using boost::math::hypot; 111 112 // The following is needed for std::complex<Real>: 113 Real fabs(const Real& a) { return abs EIGEN_NOT_A_MACRO (a); } 114 Real fmax(const Real& a, const Real& b) { using std::max; return max(a,b); } 115 116 // some specialization for the unit tests: 117 inline bool test_isMuchSmallerThan(const Real& a, const Real& b) { 118 return internal::isMuchSmallerThan(a, b, test_precision<Real>()); 119 } 120 121 inline bool test_isApprox(const Real& a, const Real& b) { 122 return internal::isApprox(a, b, test_precision<Real>()); 123 } 124 125 inline bool test_isApproxOrLessThan(const Real& a, const Real& b) { 126 return internal::isApproxOrLessThan(a, b, test_precision<Real>()); 127 } 128 129 Real get_test_precision(const Real&) { 130 return test_precision<Real>(); 131 } 132 133 Real test_relative_error(const Real &a, const Real &b) { 134 using Eigen::numext::abs2; 135 return sqrt(abs2<Real>(a-b)/Eigen::numext::mini<Real>(abs2(a),abs2(b))); 136 } 137} 138} 139 140namespace Eigen { 141 142} 143 144void test_boostmultiprec() 145{ 146 typedef Matrix<Real,Dynamic,Dynamic> Mat; 147 typedef Matrix<std::complex<Real>,Dynamic,Dynamic> MatC; 148 149 std::cout << "NumTraits<Real>::epsilon() = " << NumTraits<Real>::epsilon() << std::endl; 150 std::cout << "NumTraits<Real>::dummy_precision() = " << NumTraits<Real>::dummy_precision() << std::endl; 151 std::cout << "NumTraits<Real>::lowest() = " << NumTraits<Real>::lowest() << std::endl; 152 std::cout << "NumTraits<Real>::highest() = " << NumTraits<Real>::highest() << std::endl; 153 std::cout << "NumTraits<Real>::digits10() = " << NumTraits<Real>::digits10() << std::endl; 154 155 // chekc stream output 156 { 157 Mat A(10,10); 158 A.setRandom(); 159 std::stringstream ss; 160 ss << A; 161 } 162 { 163 MatC A(10,10); 164 A.setRandom(); 165 std::stringstream ss; 166 ss << A; 167 } 168 169 for(int i = 0; i < g_repeat; i++) { 170 int s = internal::random<int>(1,EIGEN_TEST_MAX_SIZE); 171 172 CALL_SUBTEST_1( cholesky(Mat(s,s)) ); 173 174 CALL_SUBTEST_2( lu_non_invertible<Mat>() ); 175 CALL_SUBTEST_2( lu_invertible<Mat>() ); 176 CALL_SUBTEST_2( lu_non_invertible<MatC>() ); 177 CALL_SUBTEST_2( lu_invertible<MatC>() ); 178 179 CALL_SUBTEST_3( qr(Mat(internal::random<int>(1,EIGEN_TEST_MAX_SIZE),internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) ); 180 CALL_SUBTEST_3( qr_invertible<Mat>() ); 181 182 CALL_SUBTEST_4( qr<Mat>() ); 183 CALL_SUBTEST_4( cod<Mat>() ); 184 CALL_SUBTEST_4( qr_invertible<Mat>() ); 185 186 CALL_SUBTEST_5( qr<Mat>() ); 187 CALL_SUBTEST_5( qr_invertible<Mat>() ); 188 189 CALL_SUBTEST_6( selfadjointeigensolver(Mat(s,s)) ); 190 191 CALL_SUBTEST_7( eigensolver(Mat(s,s)) ); 192 193 CALL_SUBTEST_8( generalized_eigensolver_real(Mat(s,s)) ); 194 195 TEST_SET_BUT_UNUSED_VARIABLE(s) 196 } 197 198 CALL_SUBTEST_9(( jacobisvd(Mat(internal::random<int>(EIGEN_TEST_MAX_SIZE/4, EIGEN_TEST_MAX_SIZE), internal::random<int>(EIGEN_TEST_MAX_SIZE/4, EIGEN_TEST_MAX_SIZE/2))) )); 199 CALL_SUBTEST_10(( bdcsvd(Mat(internal::random<int>(EIGEN_TEST_MAX_SIZE/4, EIGEN_TEST_MAX_SIZE), internal::random<int>(EIGEN_TEST_MAX_SIZE/4, EIGEN_TEST_MAX_SIZE/2))) )); 200} 201 202