# A constructive approach

**The problem can be solved if the form of the solution is given.**

Define the two factors using a hint (that these should be *cubic* equations) in the original post

```
y1 = x^3 - p1 x + q1
y2 = x^3 - p2 x + q2
```

Build a companion matrix of the polynomial $p(x)$

```
CompanionMatrix[p_,x_]:=Module[{n,w=CoefficientList[p,x]},w=-w/Last[w];
n=Length[w]-1; SparseArray[{{i_,n}:>w[[i]],{i_,j_}/;i==j+1->1},{n,n}]]
```

- The roots of a polynomial equation $p_A(x)=0$ are given by the
eigenvalues of its companion matrix $A$.
- If $a$ is a root of $p_A(x)$
(with companion matrix $A$) and $b$ is a root of $p_B(x)$ (with
companion matrix $B$) then $a b$ is an eigenvalue of $A\otimes B$ and
$a+b$ is an eigenvalue of $A\otimes I+I\otimes B$, where $\otimes$
is the direct product of matrices.

Let us focus on the *product* case and, therefore, determine the characteristic equation of the direct product of companion matrices

```
(a1=CompanionMatrix[y1,x])//MatrixForm
(a2=CompanionMatrix[y2,x])//MatrixForm
y3=CharacteristicPolynomial[KroneckerProduct[a1,a2],x]
y4=y3 Sign[CoefficientList[y3,x]//Last]
```

The sum is treated similarly.

We find out that the resulting characteristic polynomial is

```
y4=-q1^3 q2^3 + p1 p2 q1^2 q2^2 x + (-p2^3 q1^2 - p1^3 q2^2 + 3 q1^2 q2^2) x^3 + p1 p2 q1 q2 x^4 + p1^2 p2^2 x^5 - 3 q1 q2 x^6 - 2 p1 p2 x^7 + x^9
```

**Notice the relatively simple form in the case of trinomial equations.** Now, your polynomial is

```
z = -1 - #1 + 3 #1^3 - #1^4 + #1^5 - 3 #1^6 + 2 #1^7 + #1^9 &@x
```

We demand that the two polynomials (`y3`

and `z`

) are equal for any value of `x`

```
r = SolveAlways[z == y4, {x}]
```

Several solutions are obtained. Take, for instance, the first one

```
{y1, y2} /. r[[1]]
(*{1/q2 + x/q2^(2/3) + x^3, q2 - q2^(2/3) x + x^3}*)
```

Setting `q2=1`

we obtain the OP result.

# Comment on the method

MA has a nice `RootApproximant`

function, which operates by virtue of the LLL algorithm. One may be tempted to follow this route and implement some kind of this experimental mathematics approach. In contrast, the presented solution is *fully constructive* (in fact, algebraic) and does not require arbitrary precision computations.

# Answer to the 1st challenge question

```
y1=x^5+ m x^3+n x^2+p1 x+q1;
y2=x^3+ p2 x+q2;
a1=CompanionMatrix[y1,x];
a2=CompanionMatrix[y2,x];
y3=CharacteristicPolynomial[KroneckerProduct[a1,a2],x];
y4=y3 Sign[CoefficientList[y3,x]//Last];
z=-282300416-64550400 #-34426880 #^2-14185880 #^3+8564800 #^4+4231216 #^5-972800 #^6-367820 #^7+27360 #^8+2600 #^9+1680 #^10+100 #^11-240 #^12+40 #^13+#^15&@x;
r=SolveAlways[y4==z,{x}]
({y1,y2}/.r[[1]])/.q2->1
Out[1]= {656-150 x+80 x^2-20 x^3+x^5,1+x+x^3}
```

Notice that in the characteristic equation the highest order term can
be negative, whereas I assume that in the given equation `z`

it is
always 1. Therefore, in order to use `SolveAlways`

`y3`

is multiplied
by the sign.

From the shown solutions, one can see that there is some
arbitrariness in the results. But this is, of course, expected. As
for the shown solution in radicals, one probably needs to guess the
field extension.

One can create a `Module`

as to fully automate the derivation. But is
there a pressing need?

# Answer to the 2st challenge question

Here I demonstrate that the method can be used to write a root $P$ in terms of three roots of lower-order polynomials as $P=Q+R S$. Additionally, the computation is done in a more structured form

- Define a generic polynomial with 1 as the leading coefficient

```
pX[a_,n_]:=x^n+Sum[a[i]x^i,{i,0,n-1}]
```

- Define two modules that split a root in terms of a sum or a product

```
pSum[pA_,pB_,y_]:=Module[{mA,mB,mAB,nA,nB,p},
nA=Exponent[pA,y];
nB=Exponent[pB,y];
mA=CompanionMatrix[pA,y];
mB=CompanionMatrix[pB,y];
mAB=KroneckerProduct[mA,IdentityMatrix[nB]]+KroneckerProduct[IdentityMatrix[nA],mB];
p=CharacteristicPolynomial[mAB,y];
p Sign[CoefficientList[p,y]//Last]
]
pProduct[pA_,pB_,y_]:=Module[{mA,mB,mAB,p},
mA=CompanionMatrix[pA,y];
mB=CompanionMatrix[pB,y];
mAB=KroneckerProduct[mA,mB];
p=CharacteristicPolynomial[mAB,y];
p Sign[CoefficientList[p,y]//Last]
]
```

- Construct a working example for polynomials of the degrees 3, 2, 2
respectively.

```
{nQ,nR,nS}={3,2,2};
nT=nR nS;
pP[0]=(RootReduce[Root[#^3+#+3&,1]+Root[#^2+#+5&,1]Root[#^2+#+7&,1]]//First)@x
Out[1]= 1984051825-172403780 x-281288553 x^2+14148329 x^3+17544721 x^4-310509 x^5-619703 x^6-4623 x^7+13443 x^8+244 x^9-167 x^10-3 x^11+x^12
```

- Do the
**first part**, namely, split the root of a 12th order equation
into a sum of 3rd and 4th order roots

```
(sol[1]=SolveAlways[pP[0]==pSum[pX[q,nQ],pX[t,nT],x],x])//Transpose//TableForm
rule[1]=q[2]->0;
sol[1,1]=First[sol[1]]/.rule[1];
AppendTo[sol[1,1],rule[1]];
{pX[q,nQ],pX[t,nT]}/.sol[1,1]
Out[2]= {3+x+x^3,1225-35 x-58 x^2-x^3+x^4}
```

- Do the
**second part**, split the 4th order root into a product of 2
roots of quadratic equations

```
(sol[2]=SolveAlways[(pX[t,nT]/.sol[1,1])==pProduct[pX[r,nR],pX[s,nS],x],x])//Transpose//TableForm
rule[2]=s[1]->1;
sol[2,1]=First[sol[2]]/.rule[2];
AppendTo[sol[2,1],rule[2]];
{pX[r,nR],pX[s,nS]}/.sol[2,1]
Out[3]= {5+x+x^2,7+x+x^2}
```

- The coefficients in
`rule[1]`

and `rule[2]`

have been selected as to
match the original equation. However, other choices are possible.

5I'm sure you are aware of this but I'll say it anyway: this is not an easy problem. – Daniel Lichtblau – 2016-05-09T14:14:15.037

Would this perhaps be better asked on [Math.SE]? – Mr.Wizard – 2016-07-08T06:31:44.520

Recently we worked on a Metropolis-Hasting algorithm here in permutation-98 space. Consider two polynomials in rational coefficients of degree 3. That's rational-8 space. Is it inconceivable to adapt this algorithm to search this space and converge to a solution? – Dominic – 2019-08-07T20:49:00.127

1@Dominic for Metropolis-Hastings you need to define a reasonably continuous quality function that you can try to maximize or minimize. Do you have any ideas for how to define it? – Roman – 2019-08-07T22:57:37.927

You might get some inspiration from the Fiedler factorization of companion matrices. – yarchik – 2019-08-08T07:30:25.013

2

Just to state the obvious, going in opposite direction can be done with the help of direct products of companion matrices. Here, the problem is to write a companion matrix as a direct product of companion matrices. https://math.stackexchange.com/questions/331017/enlightening-proof-that-the-algebraic-numbers-form-a-field?noredirect=1&lq=1

– yarchik – 2019-08-08T07:35:34.090@yarchik - you can go in the opposite direction using

– Bob Hanlon – 2019-08-09T03:34:43.240`RootReduce`