Is there a way to Collect[] for more than one symbol?

Oftentimes you find yourself looking for polynomials in multiple variables. Consider the following expression:

a(x - y)^3 +  b(x - y) + c(x - y) + d

as you can see this is clearly a polynomial in x-y. Is there an equivalent of Collect, that works on more complicated expressions than just a single variable? I would like to have something similar to

Collect[%, x - y]
(* --> a(x - y)^3 + (b+c)(x - y) + d *)

however. Collect can not work on x-y. Of course you could solve this first example by substituting x-y -> z, then Collect the z, and afterwards substitute back like so:

a(x - y)^3 +  b(x - y) + c(x - y) + d /. x-y->z

gives

d + b z + c z + a z^3

Then

Collect[a z^3 + b z + c z + d, z]

gives

d + (b + c) z + a z^3

Now undo the substitution by running % /. z -> x - y. This gives the desired result:

d + (b + c) (x - y) + a (x - y)^3

So this is good. For obvious polynomials, we can solve this. But what about real world examples? Would you have guessed that

d + b x + c x + a x^3 - b y - c y - 3 a x^2 y + 3 a x y^2 - a y^3

is exactly the same polynomial? How would you Collect x-y here, as you cannot do the substitution?

• I'd make the substitution expr /. x -> y + z before applying Collect[] (and possibly Simplify[] before that) myself... – J. M. will be back soon Jan 17 '12 at 22:49
• which version are you using? This works in v.8.0.4 – rcollyer Jan 17 '12 at 22:53
• I am using mathematica 7 – niklasfi Jan 17 '12 at 23:08
• I think it has "worked" since version 3. [Quotes in use because non-variable "variables" implementation is a bit of hackery. Nested quotes in use because "'non-variable' variables" is an oxymoron. Nested quotes in use...] – Daniel Lichtblau Jan 17 '12 at 23:09
• @niklasfi: That's easily solved with a module: Module[{z},Collect[expr/.x->y+z,z]/.z->x-y]. The only limitation is that expr must not contain z. Or use Collect[expr/.x->y+#,#]/.#->x-y]&@Unique[] to lift even that limitation. – celtschk Mar 7 '12 at 8:13

Well, I am more inclined to try something that leverages Mathematica "knowledge" of polynomials.

In fact, in Mathematica 7 and 8, you can collect by $x-y$. However, it only works if $(x-y)$ is explicitly apparent in the form of the argument seen by the Collect function.

So, this works:

Collect[a (x - y)^3 + b (x - y)^2 + c (x - y) + d, x - y]

But, this doesn't:

Collect[Expand[a (x - y)^3 + b (x - y)^2 + c (x - y) + d], x - y]

I would use PolynomialReduce but I don't have an automated way for doing what you need.

Nevertheless it looks promising as the following does return {a, b, c, d}:

Flatten[
PolynomialReduce[
d + c x + b x^2 + a x^3 - c y - 2 b x y - 3 a x^2 y
+ b y^2 + 3 a x y^2 - a y^3,
Table[(x - y)^i, {i, 3, 1, -1}],
{x, y}
]
]

But I can use this approach for other factorizations.

For example, consider the following expansion Expand[z (x - y)^6 + w (x - y)^4 + t (x - y)^2 + s]. Using this expansion I am able to retrieve the coefficients in terms of powers of x^2 - 2 x y + y^2:

Flatten[
PolynomialReduce[
s + t x^2 + w x^4 - 2 t x y - 4 w x^3 y + t y^2 + 6 w x^2 y^2 -
4 w x y^3 + w y^4 + x^6 z - 6 x^5 y z + 15 x^4 y^2 z -
20 x^3 y^3 z + 15 x^2 y^4 z - 6 x y^5 z + y^6 z,
Table[(x^2 - 2 x y + y^2)^i, {i, 3, 1, -1}],
{x, y}
]
]

Well, technically you can do it programmatically:

CollectMany[expr_, a_ - b_] := Block[{cZf23},
Collect[expr /. a -> b + cZf23, cZf23] /. cZf23 -> (a - b)]

In:= CollectMany[
d + b x + c x + a x^3 - b y - c y - 3 a x^2 y + 3 a x y^2 - a y^3,
x - y]

Out= d + (b + c) (x - y) + a (x - y)^3

Notice I use the symbol cZf23. You want this to be some random and HIGHLY unlikely to be in use.

In the general case you may need a slightly more complicated definition:

CollectMany[expr_, a_ - b_] := Block[{cZf23},
Collect[expr /. a -> b + cZf23, cZf23] /. cZf23 -> (a - b)]
CollectMany[expr_, a_ + b_] := Block[{cZf23},
Collect[expr /. a -> cZf23 - b, cZf23] /. cZf23 -> (a + b)]

I'm not 100% sure if CollectMany[polynomial, a + b] will trigger in the same way as CollectMany[polynomial, a - b], so you may need to do things this way.

An updated version I would use instead is:

CollectMany[expr_, a_ + b_] := Block[{z = Unique[]},
Collect[expr /. a -> z - b, z] /. z -> (a + b)]

Using a_ + b_ will match x - y properly.

• @ThiesHeidecke: The main reason why I'm hesitant to use Module is because I've run into the unfortunate bug where defined variables inside of Module are cleared but the symbols don't disappear from the Mathematica kernel. In tight loops, this manifests as a memory leak. This question on SO addresses this. – Mike Bailey Jan 18 '12 at 0:53
• Does this happen even if you Quit[] the kernel? – rm -rf Jan 18 '12 at 6:37
• @R.M: If you quit, then the symbols get unloaded. But for a continuous session, Module will slowly leak memory. This has been a known bug since version 7. I'm not sure if it still exists at this moment. – Mike Bailey Jan 18 '12 at 8:04
• @Mike You can avoid usage of Block and Unique as well as Module simply by using a formal symbol. – Alexey Popkov Jun 18 '12 at 13:09
• @Mike In your updated version usage of the symbol z is dangerous since this symbol can be used inside of expr. – Alexey Popkov Jun 18 '12 at 13:21

An easy way to do this is:

Normal[Series[d + c x + b x^2 + a x^3 - c y - 2 b x y - 3 a x^2 y
+ b y^2 + 3 a x y^2 - a y^3, {x, y, 3}]]

(*
===> d + c (x - y) + b (x - y)^2 + a (x - y)^3
*)

This deals with what more or less is a general case: f is a polynomial in g, where g is another polynomial and we want to "collect" with respect to g. Actually, I once posted similar code on the MathGroup but could not find it so I decided to write it again. (I did not have the time to write nicer code so this is pretty rough...).

The code will only work when f actually is a polynomial in g (with constant coefficients) otherwise it will simply return f itself. I have had to use HoldForm to prevent Mathematica rearranging linear terms in some cases.

Here is the code:

CollectPolynomial[f_, g_, vars_] :=
Module[{s = PolynomialReduce[f, g, vars], u, v, coeffs, z, p},
u = s[[1, 1]]; v = s[]; coeffs = {v};
While[FreeQ[v, Alternatives @@ vars] && u =!= 0,
s = PolynomialReduce[u, g, vars]; u = s[[1, 1]]; v = s[];
coeffs = Prepend[coeffs, v]];
p = Table[z^i, {i, 0, Length[coeffs] - 1}].Reverse[coeffs] /.
z -> HoldForm[g]; If[Simplify[ReleaseHold[p] - f] === 0, p, f]]

For example:

g = x^2 + y^2 - 1;f = Expand[Sum[g^i, {i, 0, 5}]] + 3;

CollectPolynomial[f, g, {x, y}]
4 + Hold[-1 + x^2 + y^2] + Hold[-1 + x^2 + y^2]^2 +
Hold[-1 + x^2 + y^2]^3 + Hold[-1 + x^2 + y^2]^4 +
Hold[-1 + x^2 + y^2]^5