boostmultiprec.cpp (5731B)
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 #ifdef EIGEN_TEST_PART_11 59 #include "simplicial_cholesky.cpp" 60 #endif 61 62 #include <Eigen/Dense> 63 64 #undef min 65 #undef max 66 #undef isnan 67 #undef isinf 68 #undef isfinite 69 #undef I 70 71 #include <boost/serialization/nvp.hpp> 72 #include <boost/multiprecision/cpp_dec_float.hpp> 73 #include <boost/multiprecision/number.hpp> 74 #include <boost/math/special_functions.hpp> 75 #include <boost/math/complex.hpp> 76 77 namespace mp = boost::multiprecision; 78 typedef mp::number<mp::cpp_dec_float<100>, mp::et_on> Real; 79 80 namespace Eigen { 81 template<> struct NumTraits<Real> : GenericNumTraits<Real> { 82 static inline Real dummy_precision() { return 1e-50; } 83 }; 84 85 template<typename T1,typename T2,typename T3,typename T4,typename T5> 86 struct NumTraits<boost::multiprecision::detail::expression<T1,T2,T3,T4,T5> > : NumTraits<Real> {}; 87 88 template<> 89 Real test_precision<Real>() { return 1e-50; } 90 91 // needed in C++93 mode where number does not support explicit cast. 92 namespace internal { 93 template<typename NewType> 94 struct cast_impl<Real,NewType> { 95 static inline NewType run(const Real& x) { 96 return x.template convert_to<NewType>(); 97 } 98 }; 99 100 template<> 101 struct cast_impl<Real,std::complex<Real> > { 102 static inline std::complex<Real> run(const Real& x) { 103 return std::complex<Real>(x); 104 } 105 }; 106 } 107 } 108 109 namespace boost { 110 namespace multiprecision { 111 // to make ADL works as expected: 112 using boost::math::isfinite; 113 using boost::math::isnan; 114 using boost::math::isinf; 115 using boost::math::copysign; 116 using boost::math::hypot; 117 118 // The following is needed for std::complex<Real>: 119 Real fabs(const Real& a) { return abs EIGEN_NOT_A_MACRO (a); } 120 Real fmax(const Real& a, const Real& b) { using std::max; return max(a,b); } 121 122 // some specialization for the unit tests: 123 inline bool test_isMuchSmallerThan(const Real& a, const Real& b) { 124 return internal::isMuchSmallerThan(a, b, test_precision<Real>()); 125 } 126 127 inline bool test_isApprox(const Real& a, const Real& b) { 128 return internal::isApprox(a, b, test_precision<Real>()); 129 } 130 131 inline bool test_isApproxOrLessThan(const Real& a, const Real& b) { 132 return internal::isApproxOrLessThan(a, b, test_precision<Real>()); 133 } 134 135 Real get_test_precision(const Real&) { 136 return test_precision<Real>(); 137 } 138 139 Real test_relative_error(const Real &a, const Real &b) { 140 using Eigen::numext::abs2; 141 return sqrt(abs2<Real>(a-b)/Eigen::numext::mini<Real>(abs2(a),abs2(b))); 142 } 143 } 144 } 145 146 namespace Eigen { 147 148 } 149 150 EIGEN_DECLARE_TEST(boostmultiprec) 151 { 152 typedef Matrix<Real,Dynamic,Dynamic> Mat; 153 typedef Matrix<std::complex<Real>,Dynamic,Dynamic> MatC; 154 155 std::cout << "NumTraits<Real>::epsilon() = " << NumTraits<Real>::epsilon() << std::endl; 156 std::cout << "NumTraits<Real>::dummy_precision() = " << NumTraits<Real>::dummy_precision() << std::endl; 157 std::cout << "NumTraits<Real>::lowest() = " << NumTraits<Real>::lowest() << std::endl; 158 std::cout << "NumTraits<Real>::highest() = " << NumTraits<Real>::highest() << std::endl; 159 std::cout << "NumTraits<Real>::digits10() = " << NumTraits<Real>::digits10() << std::endl; 160 161 // check stream output 162 { 163 Mat A(10,10); 164 A.setRandom(); 165 std::stringstream ss; 166 ss << A; 167 } 168 { 169 MatC A(10,10); 170 A.setRandom(); 171 std::stringstream ss; 172 ss << A; 173 } 174 175 for(int i = 0; i < g_repeat; i++) { 176 int s = internal::random<int>(1,EIGEN_TEST_MAX_SIZE); 177 178 CALL_SUBTEST_1( cholesky(Mat(s,s)) ); 179 180 CALL_SUBTEST_2( lu_non_invertible<Mat>() ); 181 CALL_SUBTEST_2( lu_invertible<Mat>() ); 182 CALL_SUBTEST_2( lu_non_invertible<MatC>() ); 183 CALL_SUBTEST_2( lu_invertible<MatC>() ); 184 185 CALL_SUBTEST_3( qr(Mat(internal::random<int>(1,EIGEN_TEST_MAX_SIZE),internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) ); 186 CALL_SUBTEST_3( qr_invertible<Mat>() ); 187 188 CALL_SUBTEST_4( qr<Mat>() ); 189 CALL_SUBTEST_4( cod<Mat>() ); 190 CALL_SUBTEST_4( qr_invertible<Mat>() ); 191 192 CALL_SUBTEST_5( qr<Mat>() ); 193 CALL_SUBTEST_5( qr_invertible<Mat>() ); 194 195 CALL_SUBTEST_6( selfadjointeigensolver(Mat(s,s)) ); 196 197 CALL_SUBTEST_7( eigensolver(Mat(s,s)) ); 198 199 CALL_SUBTEST_8( generalized_eigensolver_real(Mat(s,s)) ); 200 201 TEST_SET_BUT_UNUSED_VARIABLE(s) 202 } 203 204 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))) )); 205 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))) )); 206 207 CALL_SUBTEST_11(( test_simplicial_cholesky_T<Real,int,ColMajor>() )); 208 }