boostmultiprec.cpp 5.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208
  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. #include <sstream>
  10. #ifdef EIGEN_TEST_MAX_SIZE
  11. #undef EIGEN_TEST_MAX_SIZE
  12. #endif
  13. #define EIGEN_TEST_MAX_SIZE 50
  14. #ifdef EIGEN_TEST_PART_1
  15. #include "cholesky.cpp"
  16. #endif
  17. #ifdef EIGEN_TEST_PART_2
  18. #include "lu.cpp"
  19. #endif
  20. #ifdef EIGEN_TEST_PART_3
  21. #include "qr.cpp"
  22. #endif
  23. #ifdef EIGEN_TEST_PART_4
  24. #include "qr_colpivoting.cpp"
  25. #endif
  26. #ifdef EIGEN_TEST_PART_5
  27. #include "qr_fullpivoting.cpp"
  28. #endif
  29. #ifdef EIGEN_TEST_PART_6
  30. #include "eigensolver_selfadjoint.cpp"
  31. #endif
  32. #ifdef EIGEN_TEST_PART_7
  33. #include "eigensolver_generic.cpp"
  34. #endif
  35. #ifdef EIGEN_TEST_PART_8
  36. #include "eigensolver_generalized_real.cpp"
  37. #endif
  38. #ifdef EIGEN_TEST_PART_9
  39. #include "jacobisvd.cpp"
  40. #endif
  41. #ifdef EIGEN_TEST_PART_10
  42. #include "bdcsvd.cpp"
  43. #endif
  44. #ifdef EIGEN_TEST_PART_11
  45. #include "simplicial_cholesky.cpp"
  46. #endif
  47. #include <Eigen/Dense>
  48. #undef min
  49. #undef max
  50. #undef isnan
  51. #undef isinf
  52. #undef isfinite
  53. #undef I
  54. #include <boost/serialization/nvp.hpp>
  55. #include <boost/multiprecision/cpp_dec_float.hpp>
  56. #include <boost/multiprecision/number.hpp>
  57. #include <boost/math/special_functions.hpp>
  58. #include <boost/math/complex.hpp>
  59. namespace mp = boost::multiprecision;
  60. typedef mp::number<mp::cpp_dec_float<100>, mp::et_on> Real;
  61. namespace Eigen {
  62. template<> struct NumTraits<Real> : GenericNumTraits<Real> {
  63. static inline Real dummy_precision() { return 1e-50; }
  64. };
  65. template<typename T1,typename T2,typename T3,typename T4,typename T5>
  66. struct NumTraits<boost::multiprecision::detail::expression<T1,T2,T3,T4,T5> > : NumTraits<Real> {};
  67. template<>
  68. Real test_precision<Real>() { return 1e-50; }
  69. // needed in C++93 mode where number does not support explicit cast.
  70. namespace internal {
  71. template<typename NewType>
  72. struct cast_impl<Real,NewType> {
  73. static inline NewType run(const Real& x) {
  74. return x.template convert_to<NewType>();
  75. }
  76. };
  77. template<>
  78. struct cast_impl<Real,std::complex<Real> > {
  79. static inline std::complex<Real> run(const Real& x) {
  80. return std::complex<Real>(x);
  81. }
  82. };
  83. }
  84. }
  85. namespace boost {
  86. namespace multiprecision {
  87. // to make ADL works as expected:
  88. using boost::math::isfinite;
  89. using boost::math::isnan;
  90. using boost::math::isinf;
  91. using boost::math::copysign;
  92. using boost::math::hypot;
  93. // The following is needed for std::complex<Real>:
  94. Real fabs(const Real& a) { return abs EIGEN_NOT_A_MACRO (a); }
  95. Real fmax(const Real& a, const Real& b) { using std::max; return max(a,b); }
  96. // some specialization for the unit tests:
  97. inline bool test_isMuchSmallerThan(const Real& a, const Real& b) {
  98. return internal::isMuchSmallerThan(a, b, test_precision<Real>());
  99. }
  100. inline bool test_isApprox(const Real& a, const Real& b) {
  101. return internal::isApprox(a, b, test_precision<Real>());
  102. }
  103. inline bool test_isApproxOrLessThan(const Real& a, const Real& b) {
  104. return internal::isApproxOrLessThan(a, b, test_precision<Real>());
  105. }
  106. Real get_test_precision(const Real&) {
  107. return test_precision<Real>();
  108. }
  109. Real test_relative_error(const Real &a, const Real &b) {
  110. using Eigen::numext::abs2;
  111. return sqrt(abs2<Real>(a-b)/Eigen::numext::mini<Real>(abs2(a),abs2(b)));
  112. }
  113. }
  114. }
  115. namespace Eigen {
  116. }
  117. EIGEN_DECLARE_TEST(boostmultiprec)
  118. {
  119. typedef Matrix<Real,Dynamic,Dynamic> Mat;
  120. typedef Matrix<std::complex<Real>,Dynamic,Dynamic> MatC;
  121. std::cout << "NumTraits<Real>::epsilon() = " << NumTraits<Real>::epsilon() << std::endl;
  122. std::cout << "NumTraits<Real>::dummy_precision() = " << NumTraits<Real>::dummy_precision() << std::endl;
  123. std::cout << "NumTraits<Real>::lowest() = " << NumTraits<Real>::lowest() << std::endl;
  124. std::cout << "NumTraits<Real>::highest() = " << NumTraits<Real>::highest() << std::endl;
  125. std::cout << "NumTraits<Real>::digits10() = " << NumTraits<Real>::digits10() << std::endl;
  126. // check stream output
  127. {
  128. Mat A(10,10);
  129. A.setRandom();
  130. std::stringstream ss;
  131. ss << A;
  132. }
  133. {
  134. MatC A(10,10);
  135. A.setRandom();
  136. std::stringstream ss;
  137. ss << A;
  138. }
  139. for(int i = 0; i < g_repeat; i++) {
  140. int s = internal::random<int>(1,EIGEN_TEST_MAX_SIZE);
  141. CALL_SUBTEST_1( cholesky(Mat(s,s)) );
  142. CALL_SUBTEST_2( lu_non_invertible<Mat>() );
  143. CALL_SUBTEST_2( lu_invertible<Mat>() );
  144. CALL_SUBTEST_2( lu_non_invertible<MatC>() );
  145. CALL_SUBTEST_2( lu_invertible<MatC>() );
  146. CALL_SUBTEST_3( qr(Mat(internal::random<int>(1,EIGEN_TEST_MAX_SIZE),internal::random<int>(1,EIGEN_TEST_MAX_SIZE))) );
  147. CALL_SUBTEST_3( qr_invertible<Mat>() );
  148. CALL_SUBTEST_4( qr<Mat>() );
  149. CALL_SUBTEST_4( cod<Mat>() );
  150. CALL_SUBTEST_4( qr_invertible<Mat>() );
  151. CALL_SUBTEST_5( qr<Mat>() );
  152. CALL_SUBTEST_5( qr_invertible<Mat>() );
  153. CALL_SUBTEST_6( selfadjointeigensolver(Mat(s,s)) );
  154. CALL_SUBTEST_7( eigensolver(Mat(s,s)) );
  155. CALL_SUBTEST_8( generalized_eigensolver_real(Mat(s,s)) );
  156. TEST_SET_BUT_UNUSED_VARIABLE(s)
  157. }
  158. 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))) ));
  159. 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))) ));
  160. CALL_SUBTEST_11(( test_simplicial_cholesky_T<Real,int,ColMajor>() ));
  161. }