Efficient way of finding polynomial relation

Question

I have a list p of polynomials of k, and a list of "parameters" x,y,z which are also polynomials of k.

I know that p can be written as a polynomial function of x,y,z (possibly of only a subset of the three) with integer coefficients. As an example:

p = {-2 (-1 + k - 2 k^2 + k^3), -2 (-8 - 7 k + k^3), 1 + 15 k^2 - 12 k^3};
x = {k, k + 2, 2 k};
y = {k^2 + 1, k^2 - 1, 3 k^2};
z = {1,2,3};

has a solution given by

p == 3 x^2 - 2 y x + y + 1

I need to find the relation for various p. I've been able to find some by trial and error, but the problem becomes more and more complicated as the degree of p increases. Since the lists have fixed length 3, I cannot just use Solve as I'll have an under-determined system.

As I expect the solution to be unique for each p, I've tried to construct all combinations x^a y^b z^c with degree smaller than that of p and doing a giant loop over all possible linear combinations

p == a[0] + a[1] x + a[2] y +  ... + a[nmax] x^a y^b z^c

for all coefficients a[i] in some list until I get True, but this brute force way is obviously extremely time consuming, and scales very badly with the degree.

Is there an efficient way to solve that kind of problems with Mathematica? Or at least make the trial and error process simpler.

Answer 1

You can use SolveAlways. First, I use the following ansatz:

ansatz=Sum[a[i, j, l] x^i y^j z^l, {i, 0, 2}, {j, 0, 2}, {l, 0, 1}];

where I cap the degrees of the x and y variables to 2 and the z variable to 1. Higher degrees will yield multiple results. Then, I define the following rules:

Clear[x, y, z];
rules = {
    x -> {k, k + 2, 2 k}, 
    y -> {k^2 + 1, k^2 - 1, 3 k^2},
    z -> {1, 2, 3}
};

Using SolveAlways:

sol = SolveAlways[p == ansatz /. rules, k]

{{a[0, 0, 0] -> 1, a[0, 0, 1] -> 0, a[0, 1, 0] -> 1, a[0, 1, 1] -> 0, a[0, 2, 0] -> 0, a[0, 2, 1] -> 0, a[1, 0, 0] -> 0, a[1, 0, 1] -> 0, a[1, 1, 0] -> -2, a[1, 1, 1] -> 0, a[1, 2, 0] -> 0, a[1, 2, 1] -> 0, a[2, 0, 0] -> 3, a[2, 0, 1] -> 0, a[2, 1, 0] -> 0, a[2, 1, 1] -> 0, a[2, 2, 0] -> 0, a[2, 2, 1] -> 0}}

This yields:

ansatz /. First @ sol

1 + 3 x^2 + y - 2 x y

Answer 2

Eliminate often is used for calculations like this. Define

pl = {-2 (-1 + k - 2 k^2 + k^3), -2 (-8 - 7 k + k^3), 1 + 15 k^2 - 12 k^3};
xl = {k, k + 2, 2 k};
yl = {k^2 + 1, k^2 - 1, 3 k^2};

Then, MapThread applies Eliminate to each case,

MapThread[First@Eliminate[{p == #1, x == #2, y == #3}, k] &, {pl, xl, yl}]
(* {p == -2 + 4 y - 2 x y, p == -8 + 12 x + 4 y - 2 x y, p == 1 + 5 y - 2 x y} *)

I see no general way to convert the integers to unique functions of z. In fact, even the expressions obtained for p are not necessarily unique. Consider,

MapThread[Eliminate[{p == #1, x == #2, y == #3}, {k, y}] &, {pl, xl, yl}]
(* {2 - 2 x + 4 x^2 - 2 x^3 == p, 4 - 10 x + 12 x^2 - 2 x^3 == p, 
    4 + 15 x^2 - 6 x^3 == 4 p} *)

The final case does not satisfy the integer requirement, but the first two do.