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 (at) 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 71 namespace mp = boost::multiprecision; 72 typedef mp::number<mp::cpp_dec_float<100>, mp::et_on> Real; 73 74 namespace 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 103 namespace boost { 104 namespace 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 140 namespace Eigen { 141 142 } 143 144 void 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
