Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

I've been using NSolve to find roots and now I've been wondering what underlying algorithm Mathematica uses. It's probabaly similiar to FindRoot, so is it the Newton-Raphson method?
In class we've been told that it probably is a mix of the Newton-Raphson and the bisection methods.

share|improve this question
it must be documented somewhere here specifically here and most probably there… – chris May 1 '14 at 10:01
For sparse linear systems, Solve and NSolve use several efficient numerical methods, mostly based on Gauss factoring with Markowitz products (approximately 250 pages of code). For systems of algebraic equations, NSolve computes a numerical Gröbner basis using an efficient monomial ordering, then uses eigensystem methods to extract numerical roots. See also Numerical Root Finding and Numerical Equation Solving. – Artes May 1 '14 at 11:05
For univariate inputs NSolve calls the function NRoots. That in turn has a Method option. From the documentation for NRoots: Possible settings for the Method option include: "Aberth", "CompanionMatrix", and "JenkinsTraub". The methods mentioned in the post are good for finding an individual root (not necessarily of a polynomial) but not really meant for finding all polynomial roots. And bisection is only for real roots. – Daniel Lichtblau May 1 '14 at 14:21

In at least the univariate case, the guess of "Newton-Raphson" is not too far off. Daniel mentions the three possible methods supported by NSolve[]/NRoots[] in his comment.

"CompanionMatrix" is likely done by forming the Frobenius companion matrix from the polynomial's coefficients, and then performing the usual QR algorithm to extract the eigenvalues, so that is the odd one out. (There are more efficient structure-exploiting methods than vanilla QR (e.g. this), and other possible companion matrices (e.g. this), but it doesn't seem to me that these specializations are exploited.)

The other two are basically specially-adapted Newton-Raphson methods, with somewhat better convergence properties. Jenkins-Traub, which was the default algorithm in older versions of Mathematica (I don't know now), can be interpreted as either a Newton-Raphson method applied to a specially constructed rational function, or a modified Rayleigh quotient iteration of the Frobenius companion matrix (making it a direct descendant of the classical Bernoulli method). A full description is rather involved, so I'll just link to these two papers.

Aberth's method, on the other hand, is one of a family of "simultaneous iteration" methods that are based on treating the classical Vieta formulae as a set of simultaneous nonlinear equations relating the roots and coefficients, and applying Newton-Raphson to that nonlinear system. The classical reference is this paper. (I demonstrated the derivation of the related but simpler (Weierstrass-)Durand-Kerner iteration in this math.SE answer.)

A rather comprehensive survey of all these methods is McNamee's two books.

share|improve this answer
As for the multivariate case, see this paper by Auzinger and Stetter, and this paper by a certain Mr. Daniel Lichtblau. Gröbner bases are, of course, discussed nicely in Cox/Little/O'Shea. – J. M. Jul 5 '15 at 9:42

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.