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

Question

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?

Answer 1

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}
    ]
]

Answer 2

Well, technically you can do it programmatically:

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

Using your example:

In[134]:= 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[134]= 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.

Answer 3

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
*) 

Answer 4

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[[2]]; coeffs = {v}; 
  While[FreeQ[v, Alternatives @@ vars] && u =!= 0, 
   s = PolynomialReduce[u, g, vars]; u = s[[1, 1]]; v = s[[2]]; 
   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