# Create a polynomial of a given degree

In Mathematica, how can I create a polynomial function in given variables of a given degree with unknown coefficents?

That is, I am looking for a function Poly[vars, degree] that generates, for example, if I evaluate

Poly[{x, y, z}, 3]


I should get the polynomial

a + a1*x + a2*y + a3*z + a12*xy + a13*xz + a23*yz + a11*x^2 + a22*y^2 + a33*z^2
+ a112*x^2*y + a113*x^2*z + a122*x*y^2 + a133*x*z^2 + a223*y^2*z + a233*y*z^2
+ a111*x^3 + a222*y^3 + a333*z^3


where the variables with a are unknown coefficients? Is there a simple way to do this?

• Fyi, You are missing a term. xyz there should be 20 terms. Commented Jul 16, 2016 at 4:43
• Commented Jul 16, 2016 at 18:33

polynomial[vars_List, n_Integer, coeff_] :=
#.Array[coeff, Length@#] &@ DeleteDuplicates[Times @@@ Tuples[Prepend[vars, 1], n]]

Clear[a]
polynomial[{x, y, z}, 3, a]
(* a[1] + x a[2] + y a[3] + z a[4] + x^2 a[5] + x y a[6] + x z a[7]
+ y^2 a[8] + y z a[9] + z^2 a[10] + x^3 a[11] + x^2 y a[12]
+ x^2 z a[13] + x y^2 a[14] + x y z a[15] + x z^2 a[16]
+ y^3 a[17] + y^2 z a[18] + y z^2 a[19] + z^3 a[20] *)

• I like this solution best. How would I modify it to actually have coefficients in the style "a123", as in Nasser's solution? Commented Jul 17, 2016 at 17:04
• @Posch79. Not in a convenient way. A combination of Sumit's answer and Nasser's answer would be the way to do it. Commented Jul 19, 2016 at 15:53

for a one liner

poly[{x_,y_,z_},n_,a_]:= Sum[a[i, j, k] x^i y^j z^k,
{i, 0, n}, {j, 0, n - i}, {k, 0, n - i - j}]

• The sum of the iterators should be limited to n: {i, 0, n}, {j, 0, n - i}, {k, 0, n - i - j} Commented Jul 16, 2016 at 14:03
• Thanks @BobHanlon for pointing that. Commented Jul 16, 2016 at 14:18

another way

poly[vars_List, a_, order_] := Module[{n = Length@vars, idx, z},
idx = Cases[Tuples[Range[0, order], n], x_ /; Plus @@ x <= order];
z = Times @@@ (vars^# & /@ idx);
z.((Subscript[a, Row[#]]) & /@ idx)
]

poly[{x, y, z}, a, 3]  (*a is used for coefficient*)


poly[{x, y, z}, a, 2]


poly[{x, y}, a, 2]


poly[{x}, a, 4]


poly[{x, y, z, w}, a, 5]


To make M display the coeffs first, use ParameterVariables :> {a} with TraditionalForm (this is for display only)

TraditionalForm[poly[{x, y, z}, a, 3], ParameterVariables :> {a}]


ps. If you do not like to use subscripts, you can use this instead:

poly[vars_List, order_] := Module[{n = Length@vars, idx, z},
idx = Cases[Tuples[Range[0, order], n], x_ /; Plus @@ x <= order];
z = Times @@@ (vars^# & /@ idx);
z.((ToExpression["a" <> ToString[Row[#]]]) & /@ idx)
]


and now poly[{x, y, z}, 3] gives

None of the answers thus far used one of my favorite Mathematica functions. Thus,

With[{vars = {x, y, z}, deg = 3},
Sum[With[{fs = FrobeniusSolve[ConstantArray[1, Length[vars]], k]},
Inner[#2^#1 &, fs, vars, Times].(C @@@ fs)], {k, 0, deg}]]
C[0, 0, 0] + z C[0, 0, 1] + z^2 C[0, 0, 2] + z^3 C[0, 0, 3] + y C[0, 1, 0] +
y z C[0, 1, 1] + y z^2 C[0, 1, 2] + y^2 C[0, 2, 0] + y^2 z C[0, 2, 1] + y^3 C[0, 3, 0] +
x C[1, 0, 0] + x z C[1, 0, 1] + x z^2 C[1, 0, 2] + x y C[1, 1, 0] + x y z C[1, 1, 1] +
x y^2 C[1, 2, 0] + x^2 C[2, 0, 0] + x^2 z C[2, 0, 1] + x^2 y C[2, 1, 0] + x^3 C[3, 0, 0]


If subscripts are preferred:

Instead of using ((Subscript[a, Row[#]])&/@idx) or ((ToExpression["a"<>ToString[Row[#]]])&/@idx) or Array[coeff, Length@#] or a[i, j, k] it is possible to simply use

Unique["a"]


to get nice simple unique coefficient names.

You may also use FromCoefficientRules and Indexed.

ClearAll[poly];
poly[coeff_Symbol, vars_?VectorQ, order_Integer?Positive] :=
FromCoefficientRules[
Flatten@MapIndexed[#2 - 1 -> Indexed[coeff, #2] &,
ConstantArray[0, ConstantArray[order + 1, Length@vars]], {Length@vars}],
vars]


poly builds a list of CoefficientRules using Indexed of the undefined symbol coeff to represent the coefficients. It has to use order + 1 to add space for the powers of zero since Mathematica starts counting indices at 1.

The main benefit of this method is the you can control all the coefficients through one variable and the polynomial updates automatically.

Picking a small polynomial as an example.

ClearAll[a, x, y];
p = poly[a, {x, y}, 2]


$x^2 y^2 a_{3 3}+x^2 y a_{3 2}+x^2 a_{3 1}+x y^2 a_{2 3}+x y a_{2 2}+x a_{2 1}+y^2 a_{1 3}+y a_{1 2}+a_{1 1}$

Because Indexed was used we can now set a to a $3 \times 3$ matrix (constants plus the two variables) and p will update.

SeedRandom[123];
a = RandomInteger[5, {3, 3}];
p


$3 x^2 y^2+2 x^2 y+5 x y+2 x+y^2+5 y+3$

You can change any coefficient in a and p will update. Here, the constant is changed.

a[[1, 1]] = -2;
p


$3 x^2 y^2+2 x^2 y+5 x y+2 x+y^2+5 y-2$

You can also clear the coefficients by clearing a and p will update.

a =.;
p


$x^2 y^2 a_{3 3}+x^2 y a_{3 2}+x^2 a_{3 1}+x y^2 a_{2 3}+x y a_{2 2}+x a_{2 1}+y^2 a_{1 3}+y a_{1 2}+a_{1 1}$

Hope this helps.

• This I guess demonstrates the tiny bit of friction between the choices of zero-indexing and one-indexing. :) Commented Jul 23, 2016 at 14:25