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. – Nasser Jul 16 '16 at 4:43
• – Artes Jul 16 '16 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? – Posch79 Jul 17 '16 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. – march Jul 19 '16 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} – Bob Hanlon Jul 16 '16 at 14:03
• Thanks @BobHanlon for pointing that. – Sumit Jul 16 '16 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. :) – J. M. is in limbo Jul 23 '16 at 14:25