Hi Oleg, do you mean a backport to 3.3?
gael On Wed, Mar 20, 2019 at 3:51 PM Oleg Shirokobrod <[email protected]> wrote: > Hi Gael, > > You finally modified PolynomialSolver.h code for finding roots of > polynomial with complex coefficients and wrote unit test. This enhancement > is very important for signal processing application. This modification of > PolynomialSolver.h and Companion.h are rather trivial and localized. So it > is safe to modify existing code. Would it be possible to include them to > the next release? > > Best regards, > > Oleg > > > On Wed, Nov 23, 2016 at 5:09 PM Gael Guennebaud <[email protected]> > wrote: > >> Done: >> >> https://bitbucket.org/eigen/eigen/commits/3842186ae759/ >> https://bitbucket.org/eigen/eigen/commits/eff1bfef71a2/ >> https://bitbucket.org/eigen/eigen/commits/df6c723906fc/ >> >> gael >> >> On Wed, Nov 23, 2016 at 3:39 PM, Gael Guennebaud < >> [email protected]> wrote: >> >>> Sorry for late reply, but thanks for the patches. >>> >>> Regarding the eigen-solver choice issue, it's possible to automatically >>> use the right one based on NumTraits<Scalar>::IsComplex and conditional<>. >>> >>> We would also need to extended the unit test to cover complexes. >>> >>> gael >>> >>> On Tue, Nov 8, 2016 at 3:41 PM, Oleg Shirokobrod < >>> [email protected]> wrote: >>> >>>> Dear list, >>>> >>>> Polynomial solver currently finds roots only for real polynomial >>>> coefficients. Having Eigen::ComplexEigenSolver lets Polynomial module to >>>> compute roots for complex polynomial coefficients. I have made necessary >>>> modifications. I have attached corresponding patch files. Unfortunately >>>> module cannot deduce type of the solver from partial specialization of >>>> corresponding matrix type, such that >>>> >>>> template<typename RealScalar> class EigenSolver {}; >>>> >>>> template<typenameRealScalar> class >>>> EigenSolver<Eigen::Matrix<RealScalar,_Deg,_Deg> > {// EigenSolver }; >>>> >>>> template<typename RealScalar> class >>>> EigenSolver<Eigen::Matrix<std::complex<RealScalar>,_Deg,_Deg > > { >>>> //ComplexEigenSolver}; >>>> >>>> So I have to add additional template parameter EigenSolverType with >>>> default value for real coefficients: >>>> >>>> template< typename _Scalar, int _Deg , typename EigenSolverType = >>>> EigenSolver<Matrix<_Scalar,_Deg,_Deg> > > >>>> class PolynomialSolver : public PolynomialSolverBase<_Scalar,_Deg>{}; >>>> >>>> I have to replace in number of places in the file Companion.h (mainly >>>> in functions balanced() and balance()) Scalar with RealScalar, where >>>> variables are really real. >>>> >>>> With this code I have run test against Matlab and it gave similar >>>> results. >>>> >>>> Test >>>> VectorXcd roots = VectorXcd::Random(4); >>>> VectorXcd polynomialCoefficients; >>>> roots_to_monicPolynomial(roots, polynomialCoefficients); >>>> >>>> cout << "roots : " << endl; >>>> cout << setprecision(14) << roots << endl; >>>> cout << "polynomialCoefficients : " << endl; >>>> cout << setprecision(14) << polynomialCoefficients << endl; >>>> >>>> PolynomialSolver<std::complex<double>, Dynamic, >>>> ComplexEigenSolver<MatrixXcd> > psolve(polynomialCoefficients); >>>> >>>> VectorXcd computedRoots = psolve.roots(); >>>> >>>> cout << "computedRoots : " << endl; >>>> cout << setprecision(14) << computedRoots << endl; >>>> >>>> for(int index = 0; index < computedRoots.size(); ++index) >>>> { >>>> cout << setprecision(14) << "polynom at computedRoots * 1.0 e-11: >>>> " << poly_eval(polynomialCoefficients, computedRoots(index))*1.0e11 << >>>> endl; >>>> } >>>> >>>> Output: >>>> roots : >>>> (0.12717062898648,-0.99749748222297) >>>> (0.61748100222785,-0.61339152195807) >>>> (-0.04025391399884,0.17001861629078) >>>> (0.79192480239265,-0.29941709646901) >>>> >>>> polynomialCoefficients : >>>> (0.091649023983597,-0.091441416775918) >>>> (0.24020598856634,0.37401934653925) >>>> (-0.16301627948124,-1.8544616197629) >>>> (-1.4963225196081,1.7402874843593) >>>> (1,0) >>>> >>>> computedRoots : >>>> (-0.04025391399884,0.17001861629078) >>>> (0.79192480239265,-0.29941709646901) >>>> (0.61748100222785,-0.61339152195807) >>>> (0.12717062898648,-0.99749748222297) >>>> >>>> polynom at computedRoots * 1.0 e-11: >>>> (8.3266726846887e-006,-1.2490009027033e-005) >>>> polynom at computedRoots * 1.0 e-11: >>>> (4.7184478546569e-005,-1.6653345369377e-005) >>>> polynom at computedRoots * 1.0 e-11: >>>> (2.2204460492503e-005,-1.3877787807814e-005) >>>> polynom at computedRoots * 1.0 e-11: >>>> (1.5163809063462e-005,-2.7286889655471e-005) >>>> >>>> Matlab for the same coefficients gives following roots (remember that >>>> matlab array of coefficients and Eigen one are reversed to each other) >>>> >>>> computedRoots : >>>> 0.127170628986480 - 0.997497482222969i >>>> 0.617481002227849 - 0.613391521958071i >>>> 0.791924802392650 - 0.299417096469009i >>>> -0.040253913998840 + 0.170018616290780i >>>> >>>> Best regards, >>>> >>>> Oleg Shirokobrod >>>> >>>> >>>> >>> >>
