# Factor a polynomial over the reals

I can't seem to find a function in Mathematica that factors a polynomial over the reals. Obviously, Factor doesn't work since it works over the integers, over $\mathbb{Z}_p$, or some algebraic fields extensions of the type $\mathbb{Q}(\alpha_1, \ldots, \alpha_n)$ where $\alpha_1, \ldots, \alpha_n$ are algebraic numbers. My question is general, although it arouse from my need to factor over $\mathbb{R}$ the polynomial $$f(\lambda) = \lambda^4 + 17\lambda^2+12.$$ Of course, in this particular case I can find the roots and factor it myself.

One way is to find the roots, separate into real and complex, further separate the complex ones into conjugate pairs, then reform as a factorization (taking into account the leading coefficient). I'll illustrate with the given example.

poly = x^4 + 17 x^2 + 12;
roots = x /. NSolve[poly];
rroots = Select[rts, FreeQ[#, Complex] &]
croots = Complement[roots, rroots];
topquad = Select[croots, Im[#] > 0 &]
mult = Coefficient[poly, x, Exponent[poly, x]]
rfax = x - rroots

(* Out[54]= {}

Out[56]= {0. + 0.859018 I, 0. + 4.03263 I}

Out[57]= 1

Out[58]= {}

Out[59]= {0.737913 + x^2, 16.2621 + x^2} *)


Putting together the result:

Apply[Times, Join[{mult}, rfax, cfax]]

(* Out[61]= (0.737913 + x^2) (16.2621 + x^2) *)

• After checking for irreducibility (IrreduciblePolynomialQ[]) it would be "easier" to use NumberFieldSignature[] to get the counts of real roots and complex pairs of roots. Sep 5, 2016 at 18:22
• @EricTowers If working in the setting of approximate numbers and speed is an issue then I would advocate against that unless degree is large enough to make the numerics fail. Even then it is not clear what is the gain to knowing the counts-- one needs actual root values. Sep 5, 2016 at 19:23
• Thank you, this is the nicest numerical solution.Of course it would work also if Mathematica can find the complex roots exactly. Sep 6, 2016 at 12:55
• You can replace NSolve with Solve to get exact roots. To obtain real ones, change the Select predicate to FreeQ[N[#], Complex] & and similar for separating complex roots with positive vs negative imaginary components. Sep 6, 2016 at 14:28

My knowledge in this area is quite lacking, so there are probably better ways. The following function will force-factor anything, and will return Root objects if necessary.

factor[poly_, x_] :=
Module[{n, nreal},
n = Exponent[poly, x];
nreal = CountRoots[poly, x];
Times @@ Join[
Table[(x - ToRadicals@Root[poly & /. x -> #, i]), {i, nreal}],
Table[
With[{r1 = ToRadicals@Root[poly & /. x -> #, i],
r2 = ToRadicals@Root[poly & /. x -> #, i + 1]},
(x^2 - Expand[(r1 + r2)] x + Expand[r1 r2])
],
{i, nreal + 1, n, 2}
]
]
]


This is relying on the ordering used by Root, described under "Details". Real roots come first, then complex conjugate pairs.

Example:

factor[x^4 + 17 x^2 + 12, x]
(* (17/2 - Sqrt[241]/2 + x^2) (17/2 + Sqrt[241]/2 + x^2) *)


Any coefficient in the result will be real, regardless of whether it includes terms that look imaginary in its expression. Example:

factor[x^4 + 1, x]
(* (1 - ((-1)^(1/4) - (-1)^(3/4)) x + x^2) (1 - (-(-1)^(1/4) + (-1)^(3/4)) x + x^2) *)


(-1)^(1/4) has an imaginary part, but (-1)^(1/4) + (-1)^(3/4), which is the actual coefficient of x, doesn't.

Expect not very useful results for cases like this:

factor[x^5 + 2 x^4 + x^3 + x^2 + x + 1, x]
(* (x - Root[1 + #1 + #1^2 + #1^3 + 2 #1^4 + #1^5 &, 1]) (x^2 +
Root[1 + #1 + #1^2 + #1^3 + 2 #1^4 + #1^5 &, 2] Root[
1 + #1 + #1^2 + #1^3 + 2 #1^4 + #1^5 &, 3] -
x (Root[1 + #1 + #1^2 + #1^3 + 2 #1^4 + #1^5 &, 2] +
Root[1 + #1 + #1^2 + #1^3 + 2 #1^4 + #1^5 &, 3])) (x^2 +
Root[1 + #1 + #1^2 + #1^3 + 2 #1^4 + #1^5 &, 4] Root[
1 + #1 + #1^2 + #1^3 + 2 #1^4 + #1^5 &, 5] -
x (Root[1 + #1 + #1^2 + #1^3 + 2 #1^4 + #1^5 &, 4] +
Root[1 + #1 + #1^2 + #1^3 + 2 #1^4 + #1^5 &, 5])) *)


This is why yarchik's solution is better.

• Thank you, Szabolcs. As expected, for the polynomial $f_1(\lambda) = 180+1027 \lambda ^2+666 \lambda ^4+27 \lambda ^6$ I got bad results :) By the way, your module works for monic polynomials. Sep 5, 2016 at 11:25
• Why not change Expand to ComplexExpand? Then factor[x^4 + 1, x] yields (1 - Sqrt[2] x + x^2) (1 + Sqrt[2] x + x^2) Apr 17, 2020 at 16:10
Factor[x^4 + 17 x^2 + 12, Extension -> Sqrt[241]]
(*-(1/4) (-17 + Sqrt[241] - 2 x^2) (17 + Sqrt[241] + 2 x^2)*)


The answer to another example can be obtained as follows:

f = 180 + 1027 x^2 + 666 x^4 + 27 x^6
Factor[f, Extension -> { (-293657 + 20 I Sqrt[4568243])^(1/3)}] /. (-293657 + 20 I Sqrt[4568243]) -> a // Simplify


Finally a small remark on mathematical side. Sometimes algebraic equations of higher orders (>4) can be solved in radicals. Otherwise solutions for $n=5,6,7$ are known in terms of more general functions, see for instance the wiki article on the septic equation. Very simple solutions are known for trinomial equations $x^n+x-q=0$ in terms of hypergeometric functions ($n\in\mathbb{N}$). As far as I know, these solutions and more general solutions are not directly implemented in MA. Your two examples are of course solvable because of the substitution $\lambda^2=x$.

• Thank you, but this answers only the particular question. If I have a polynomial of greater degree, this solution wouldn't work. Sep 5, 2016 at 10:02
• @Veliko Please, show your polynomial of greater degree then. There is no guarantee i will work in general, but let us try. Sep 5, 2016 at 10:03
• @Veliko What is there are no roots that can be expressed as radicals? Would you still want it factored in terms of Root expressions? Sep 5, 2016 at 10:08
• @yarchik I expect that my polynomials can be expressed as radicals. Otherwise I don't expect a Mathematica solution:) For example, other polynomials are : $f_1(\lambda) = 180 + 1027 \lambda^2 + 666 \lambda^4 + 27 \lambda^6$ and $f_2(\lambda) = 357696+5872364 \lambda ^2+13665873 \lambda ^4+8310318 \lambda ^6+573049 \lambda ^8$. Sep 5, 2016 at 11:20
• @Veliko for your first example see my updated answer. Sep 5, 2016 at 11:42