devtools/eigen
1
0
Fork 0
mirror of https://gitlab.com/libeigen/eigen.git synced 2026-04-10 11:34:33 +08:00
Code Issues Packages Projects Releases Wiki Activity

Various fixes in polynomial solver and its unit tests:

Browse Source 
- cleanup noise in imaginary part of real roots
 - take into account the magnitude of the derivative to check roots.
 - use <= instead of < at appropriate places
(grafted from 450dc97c6b
)
This commit is contained in:
Gael Guennebaud
2018-12-09 22:54:39 +01:00
parent 8fb28db12d
commit 7c42084503
3 changed files with 46 additions and 8 deletions
Download Patch File Download Diff File Expand all files Collapse all files

21
unsupported/test/polynomialsolver.cpp
View File

@@ -26,6 +26,16 @@ struct increment_if_fixed_size
}
}
template<typename PolynomialType>
PolynomialType polyder(const PolynomialType& p)
{
typedef typename PolynomialType::Scalar Scalar;
PolynomialType res(p.size());
for(Index i=1; i<p.size(); ++i)
res[i-1] = p[i]*Scalar(i);
res[p.size()-1] = 0.;
return res;
}
template<int Deg, typename POLYNOMIAL, typename SOLVER>
bool aux_evalSolver( const POLYNOMIAL& pols, SOLVER& psolve )
@@ -44,10 +54,17 @@ bool aux_evalSolver( const POLYNOMIAL& pols, SOLVER& psolve )
psolve.compute( pols );
const RootsType& roots( psolve.roots() );
EvalRootsType evr( deg );
POLYNOMIAL pols_der = polyder(pols);
EvalRootsType der( deg );
for( int i=0; i<roots.size(); ++i ){
evr[i] = std::abs( poly_eval( pols, roots[i] ) ); }
evr[i] = std::abs( poly_eval( pols, roots[i] ) );
der[i] = numext::maxi(RealScalar(1.), std::abs( poly_eval( pols_der, roots[i] ) ));
}
bool evalToZero = evr.isZero( test_precision<Scalar>() );
// we need to divide by the magnitude of the derivative because
// with a high derivative is very small error in the value of the root
// yiels a very large error in the polynomial evaluation.
bool evalToZero = (evr.cwiseQuotient(der)).isZero( test_precision<Scalar>() );
if( !evalToZero )
{
cerr << "WRONG root: " << endl;