Removing certain terms from (matrix of) n-variable polynomials

Question

Inspired by e.g. this, hopefully an easy question:

I have an n x n matrix consisting of polynomials in k + m variables x1,...,xm and y1,...,yk. In each entry of the matrix, I'd like to only keep those terms constant * x1^i1 * ... * xm^im * y1^j1 * ... * yk^jk with at most one nonzero entry among the i1,...,im, and at most one nonzero entry among the j1,...,jk, and this entry can only be 0 or 1. So I just need the total degree in the x-component to be 1 at most, and the total degree in the y-component to be 1 at most. I.e. from the polynomial

4 + 8 x1 + 3 x1 x2 + 5 x1 y4 + 8 x2^2 + 6 x1 y2^3

it should extract

4 + 8 x1 + 5 x1 y4.

So the term 3 x1 x2 is discarded because there are two nonzero powers in the i's appearing, the term 8 x2^2 is discarded because the power appearing is greater than 1, and 6 x1 y2^3 is discarded for the same reason. And 5 x1 y4 is not discarded because the two powers appearing are in the separate i's and j's.

My current guess is that I should set up something using FromCoefficientRules, and then apply that to my matrix M using Map... any help or tips would be greatly appreciated.

Answer 1

For this one can just have replacement rules that force all products (including squares) within a variable class to be zero.

xvars = Array[x, 4];
yvars = Array[y, 3];
xrels = Union[Flatten[Outer[Times, xvars, xvars]]];
yrels = Union[Flatten[Outer[Times, yvars, yvars]]];
vars = Join[xvars, yvars];
rels = Join[xrels, yrels];

To obtain a nontrivial example (one that won't give all zeros at the end), we'll use random polynomials that tend to weight the monomial powers in a way that favors lower degrees.

randpoly[vars_, len_] := 
 Module[{vlen = Length[vars]}, 
  RandomInteger[{-10, 10}, len] . 
   Table[Apply[Times, 
     RandomChoice[vars, 
      RandomChoice[Reverse[Range[vlen]] -> Range[vlen]]]], {len}]]

Now create a 3 x 3 matrix with polynomials of length up to 8 (could be less, if there are repeated terms).

SeedRandom[1234];
mat = Table[randpoly[Join[xvars, yvars], 8],
  {3}, {3}]

(* {{-6 x[2]^2 x[3] + 7 x[1] x[3]^2 + 10 x[2] x[4]^2 + 10 x[4] y[1]^2 - 
   18 y[2] + 6 x[2] y[3] - 10 x[1] y[1] y[3], -10 x[1] + 
   9 x[1] y[1] - 8 x[3] x[4] y[1] - 3 x[4]^2 y[2] - 9 y[3] - 
   2 x[2] x[3] y[3] - 3 x[1] x[3] y[1] y[2]^2 y[3], 
  2 x[2] + 10 x[3]^2 - 10 x[1] y[1]^4 + 9 x[1] x[2]^2 y[1] y[2] - 
   6 x[4]^2 y[1] y[2] - 5 y[3] - 
   3 x[2] x[3] y[1] y[3]}, {10 x[1]^3 x[2] - 4 x[3]^2 x[4] - 
   8 x[3] x[4]^2 - y[1] + 5 y[2] + 3 y[3] - x[3] y[3] + 
   y[1] y[2] y[3], -2 y[1] + 4 x[2]^2 y[1]^3 + 6 x[1] x[4] y[1]^3 - 
   8 x[4] y[2] + 2 x[2] x[3] y[2]^2 + 7 x[1]^2 x[2] x[4] y[3] - 
   y[3]^2 - 4 x[1] x[4]^2 y[1] y[3]^3, 
  5 x[2] + 10 x[1] x[3] - 5 x[1]^2 x[4] + 6 x[3] x[4] - 9 x[2] y[1] - 
   8 x[3] y[1]^2 - 6 x[2] y[3] + 8 x[1] x[3] y[3]}, {-6 x[3] + 
   10 x[4] + 5 y[3] + 7 x[1] y[3] - 6 x[3] y[3] - x[1] y[1] y[3] + 
   10 x[2] x[3] y[1] y[3] + 4 x[1]^2 y[1]^2 y[3]^2, 
  4 x[1] + 6 x[3] + 10 x[1] x[3] - 3 x[1] x[3]^2 y[1] - 7 y[2] - 
   x[3]^2 x[4] y[2]^2 - 10 x[2] y[3] + 8 x[1] x[3] x[4]^3 y[3], 
  9 x[2] - 7 x[4] + 10 y[1] - 5 x[3]^2 y[1] - 9 x[2] y[1]^2 + 
   6 x[1] x[2] x[4] y[2] + 2 x[4] y[3] - 2 x[4] y[1] y[3]}} *)

We use polynomial reduction by the relations to get rid of all bad products and powers. This will leave us with linear terms and quadratics that are products of a single variable from each class.

redmat0 = PolynomialReduce[mat, rels, vars]
  [[All, All, 2]]

(* Out[180]= {{-18 y[2] + 6 x[2] y[3], -10 x[1] + 9 x[1] y[1] - 9 y[3], 
  2 x[2] - 5 y[3]}, {-y[1] + 5 y[2] + 3 y[3] - x[3] y[3], -2 y[1] - 
   8 x[4] y[2], 
  5 x[2] - 9 x[2] y[1] - 6 x[2] y[3]}, {-6 x[3] + 10 x[4] + 5 y[3] + 
   7 x[1] y[3] - 6 x[3] y[3], 4 x[1] + 6 x[3] - 7 y[2] - 10 x[2] y[3],
   9 x[2] - 7 x[4] + 10 y[1] + 2 x[4] y[3]}} *)

We use a standard tactic to keep the quadratics and remove the linear terms.

redmat = redmat0 /. Thread[vars -> t*vars] /.
  t^2 -> 1 /. t -> 0

(* Out[181]= {{6 x[2] y[3], 9 x[1] y[1], 
  0}, {-x[3] y[3], -8 x[4] y[2], -9 x[2] y[1] - 
   6 x[2] y[3]}, {7 x[1] y[3] - 6 x[3] y[3], -10 x[2] y[3], 
  2 x[4] y[3]}} *)

Answer 2

This seems a bit hackish, and it isn't very extensible if your conditions change. Let's define a function to determine if a term is "simple".

IsSimpleTerm[term_] :=
  And[
    FreeQ[term, Power],
    DuplicateFreeQ[StringTake[SymbolName /@ Variables[term], 1]]]

Now, a "simple" polynomial is one where every term is "simple".

IsSimplePoly[poly_] := AllTrue[poly, IsSimpleTerm]

I'm not sure what you mean by "keeping terms" in your matrix. If we just delete terms that don't satisfy our condition, we probably won't have a rectangular matrix anymore. But you can map our predicate over your matrix, and maybe later use Pick.

Map[IsSimplePoly, yourMatrix, {2}]

Pick[yourMatrix, Map[IsSimplePoly, yourMatrix, {2}]]

Obviously, being able to use FreeQ[Power] is just lucky. An alternative to matching the variables on the first letter might be to create lists of the variable names that are "related" and make sure you only have one from each list. Anyway, I guess I hope you don't have to extend this :) .

Answer 3

Inspired by Daniel Lichtblau's answer, how about:

xvars = {x1, x2}
yvars = {y1, y2, y3, y4}
ruleST = Join[Thread[xvars -> t xvars], Thread[yvars -> s yvars]]
ruleTerms = {t^k_ /; k > 1 :> 0, s^k_ /; k > 1 :> 0, s -> 1, t -> 1}

For the given scalar example

mat = 4 + 8 x1 + 3 x1 x2 + 5 x1 y4 + 8 x2^2 + 6 x1 y2^3

the substitution

mat /. ruleST /. ruleTerms

produces the desired answer

4 + 8 x1 + 5 x1 y4

It works without modification if mat is an array.

Alternatively, if you want to use polynomial reduction, then

rels = {t^2, s^2, t - 1, s - 1}
vars = {s, t}
PolynomialReduce[mat /. ruleST, rels, vars][[2]]

gives the same answer. For matrices use

PolynomialReduce[mat /. ruleST, rels, vars][[All,All,2]]

as suggested by Daniel Lichtblaus.