1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2010 Manuel Yguel <manuel.yguel (at) gmail.com> 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 "main.h" 11 #include <unsupported/Eigen/Polynomials> 12 #include <iostream> 13 14 using namespace std; 15 16 namespace Eigen { 17 namespace internal { 18 template<int Size> 19 struct increment_if_fixed_size 20 { 21 enum { 22 ret = (Size == Dynamic) ? Dynamic : Size+1 23 }; 24 }; 25 } 26 } 27 28 template<typename _Scalar, int _Deg> 29 void realRoots_to_monicPolynomial_test(int deg) 30 { 31 typedef internal::increment_if_fixed_size<_Deg> Dim; 32 typedef Matrix<_Scalar,Dim::ret,1> PolynomialType; 33 typedef Matrix<_Scalar,_Deg,1> EvalRootsType; 34 35 PolynomialType pols(deg+1); 36 EvalRootsType roots = EvalRootsType::Random(deg); 37 roots_to_monicPolynomial( roots, pols ); 38 39 EvalRootsType evr( deg ); 40 for( int i=0; i<roots.size(); ++i ){ 41 evr[i] = std::abs( poly_eval( pols, roots[i] ) ); } 42 43 bool evalToZero = evr.isZero( test_precision<_Scalar>() ); 44 if( !evalToZero ){ 45 cerr << evr.transpose() << endl; } 46 VERIFY( evalToZero ); 47 } 48 49 template<typename _Scalar> void realRoots_to_monicPolynomial_scalar() 50 { 51 CALL_SUBTEST_2( (realRoots_to_monicPolynomial_test<_Scalar,2>(2)) ); 52 CALL_SUBTEST_3( (realRoots_to_monicPolynomial_test<_Scalar,3>(3)) ); 53 CALL_SUBTEST_4( (realRoots_to_monicPolynomial_test<_Scalar,4>(4)) ); 54 CALL_SUBTEST_5( (realRoots_to_monicPolynomial_test<_Scalar,5>(5)) ); 55 CALL_SUBTEST_6( (realRoots_to_monicPolynomial_test<_Scalar,6>(6)) ); 56 CALL_SUBTEST_7( (realRoots_to_monicPolynomial_test<_Scalar,7>(7)) ); 57 CALL_SUBTEST_8( (realRoots_to_monicPolynomial_test<_Scalar,17>(17)) ); 58 59 CALL_SUBTEST_9( (realRoots_to_monicPolynomial_test<_Scalar,Dynamic>( 60 internal::random<int>(18,26) )) ); 61 } 62 63 64 65 66 template<typename _Scalar, int _Deg> 67 void CauchyBounds(int deg) 68 { 69 typedef internal::increment_if_fixed_size<_Deg> Dim; 70 typedef Matrix<_Scalar,Dim::ret,1> PolynomialType; 71 typedef Matrix<_Scalar,_Deg,1> EvalRootsType; 72 73 PolynomialType pols(deg+1); 74 EvalRootsType roots = EvalRootsType::Random(deg); 75 roots_to_monicPolynomial( roots, pols ); 76 _Scalar M = cauchy_max_bound( pols ); 77 _Scalar m = cauchy_min_bound( pols ); 78 _Scalar Max = roots.array().abs().maxCoeff(); 79 _Scalar min = roots.array().abs().minCoeff(); 80 bool eval = (M >= Max) && (m <= min); 81 if( !eval ) 82 { 83 cerr << "Roots: " << roots << endl; 84 cerr << "Bounds: (" << m << ", " << M << ")" << endl; 85 cerr << "Min,Max: (" << min << ", " << Max << ")" << endl; 86 } 87 VERIFY( eval ); 88 } 89 90 template<typename _Scalar> void CauchyBounds_scalar() 91 { 92 CALL_SUBTEST_2( (CauchyBounds<_Scalar,2>(2)) ); 93 CALL_SUBTEST_3( (CauchyBounds<_Scalar,3>(3)) ); 94 CALL_SUBTEST_4( (CauchyBounds<_Scalar,4>(4)) ); 95 CALL_SUBTEST_5( (CauchyBounds<_Scalar,5>(5)) ); 96 CALL_SUBTEST_6( (CauchyBounds<_Scalar,6>(6)) ); 97 CALL_SUBTEST_7( (CauchyBounds<_Scalar,7>(7)) ); 98 CALL_SUBTEST_8( (CauchyBounds<_Scalar,17>(17)) ); 99 100 CALL_SUBTEST_9( (CauchyBounds<_Scalar,Dynamic>( 101 internal::random<int>(18,26) )) ); 102 } 103 104 void test_polynomialutils() 105 { 106 for(int i = 0; i < g_repeat; i++) 107 { 108 realRoots_to_monicPolynomial_scalar<double>(); 109 realRoots_to_monicPolynomial_scalar<float>(); 110 CauchyBounds_scalar<double>(); 111 CauchyBounds_scalar<float>(); 112 } 113 } 114
