Resultant of two polynomials

Question

The resultant of two polynomials is a polynomial in their coefficients that is zero if and only if \$p\$ and \$q\$ have a common root. It is a useful tool for eliminating variables from systems of polynomial equations. For example, it is used in my answer to the code golf challenge Add up two algebraic numbers.

The resultant of two polynomials \$p(x)\$ and \$q(x)\$ is defined as the determinant of the Sylvester matrix. More precisely, if \$p(x) = a_n x^n + a_{n-1} x^{n-1} + \cdots + a_0\$ and \$q(x) = b_m x^m + b_{m-1} x^{m-1} + \cdots + b_0\$, then the Sylvester matrix is the following \$(n+m) \times (n+m)\$ matrix:

$$ \begin{bmatrix} a_n & a_{n-1} & \cdots & a_0 & 0 & \cdots & 0 \\ 0 & a_n & a_{n-1} & \cdots & a_0 & \cdots & 0 \\ \vdots & & & & & & \vdots \\ 0 & \cdots & 0 & a_n & a_{n-1} & \cdots & a_0 \\ b_m & b_{m-1} & \cdots & b_0 & 0 & \cdots & 0 \\ 0 & b_m & b_{m-1} & \cdots & b_0 & \cdots & 0 \\ \vdots & & & & & & \vdots \\ 0 & \cdots & 0 & b_m & b_{m-1} & \cdots & b_0 \end{bmatrix} $$

where the first \$m\$ rows are the coefficients of \$p(x)\$ and the last \$n\$ rows are the coefficients of \$q(x)\$. The resultant is the determinant of this matrix.

For example, the resultant of \$p(x) = x^3 + 2x^2 + 3x + 4\$ and \$q(x) = 5x^2 + 6x + 7\$ is

$$ \begin{vmatrix} 1 & 2 & 3 & 4 & 0 \\ 0 & 1 & 2 & 3 & 4 \\ 5 & 6 & 7 & 0 & 0 \\ 0 & 5 & 6 & 7 & 0 \\ 0 & 0 & 5 & 6 & 7 \end{vmatrix} = 832. $$

In particular, when \$m = n = 0\$, the Sylvester matrix is the empty matrix, whose determinant is \$1\$.

Task

Given two nonzero polynomials \$p(x)\$ and \$q(x)\$, compute their resultant.

The coefficients of the polynomials are all integers.

You may take the polynomial in any reasonable format. For example, the polynomial \$x^4-4x^3+5x^2-2x\$ may be represented as:

• a list of coefficients, in descending order: [1,-4,5,-2,0];
• a list of coefficients, in ascending order:[0,-2,5,-4,1];
• a string representation of the polynomial, with a chosen variable, e.g., x: "x^4-4*x^3+5*x^2-2*x";
• a built-in polynomial object, e.g., x^4-4*x^3+5*x^2-2*x in PARI/GP.

When you take input as a list of coefficients, the leading coefficient is guaranteed to be nonzero.

This is code-golf, so the shortest code in bytes wins.

Test cases

Here I use coefficient lists in descending order:

[1] [2] -> 1
[3,3] [1,2,1] -> 0
[1,3,3,1] [1,0,-1] -> 0
[1,2,3,4] [5,6,7] -> 832
[1,2,3,4] [4,3,2,1] -> -2000
[1,-4,5,-2,0] [1,-4,5,-2,1] -> 1
[1,-4,5,-2,0] [1,-12,60,-160,240,-192,64,0] -> 0

Answer 1

Python NumPy, 36 bytes (-19 thanks to @DanielSchepler)

lambda p,q:p.c[0]**q.o*q(p.r).prod()

Attempt This Online!

Python NumPy, 55 bytes

lambda p,q:(p.r-q.r[:,1>0]).prod()*(p**q.o*q**p.o).c[0]

Attempt This Online!

Takes two poly1d objects as input. Has floating point issues.

Without poly1d

Python NumPy, 87 bytes

lambda p,q:prod(roots(p)-c_[roots(q)])*p[0]**~-len(q)*q[0]**~-len(p)
from numpy import*

Attempt This Online!

Takes two sequences. Has floating point issues.

How?

Uses the well-known (?) formula

\$\mathrm{res}(p,q)=p_m^nq_n^m\prod(\lambda_i-\mu_j)\$

or its variant (thanks @DanielSchepler)

\$\mathrm{res}(p,q)=p_m^n \prod q(\lambda_i)=(-1)^{mn}q_n^m \prod p(\mu_j)\$

that expresses the resultant in terms of the roots \$\{\lambda_i\}\$,\$\{\mu_j\}\$ degrees \$m,n\$ and leading coefficients of the two polynomials \$p,q\$.

Answer 2

C (gcc), 176 bytes

float f(float*p,int n,float*q,int m){return n<=m?m?*p?({int i=1;for(;i<=n;++i)q[i]-=*q/ *p*p[i];}),*p*f(p,n,q+1,m-1):n?(m%2?-*q:*q)*f(p+1,n-1,q,m):0:1:(n*m%2?-1:1)*f(q,m,p,n);}

(Note that this uses the gcc extension for block expressions, in order to be able to avoid repeated return keywords and use ?: instead of if/else.)

The idea of this submission is to perform an algorithm very reminiscent of the Euclidean algorithm on polynomials. In terms of the Sylvester matrix, you could view the reduction for example as: suppose you start with the matrix \$\$\begin{bmatrix} 5 & 6 & 7 & 0 & 0 \\ 0 & 5 & 6 & 7 & 0 \\ 0 & 0 & 5 & 6 & 7 \\ 1 & 2 & 3 & 4 & 0 \\ 0 & 1 & 2 & 3 & 4\end{bmatrix}.\$\$ Then by subtracting one fifth of row 1 from row 4 and also subtracting one fifth of row 2 from row 5, you reduce to taking the determinant of: \$\$\begin{bmatrix} 5 & 6 & 7 & 0 & 0 \\ 0 & 5 & 6 & 7 & 0 \\ 0 & 0 & 5 & 6 & 7 \\ 0 & 4/5 & 8/5 & 4 & 0 \\ 0 & 0 & 4/5 & 8/5 & 4 \end{bmatrix}.\$\$ Now, by expanding by minors along the first column, this is 5 times the resultant of \$5x^2 + 6x + 7\$ and \$(4/5) x^2 + (8/5) x + 4\$.

(Do note that in the intermediate steps, it is possible that it is considering a resultant of polynomials with one of the leading coefficients being 0. In this case, it is still calculating the determinant of a generalized Sylvester matrix.)

(It's also slightly amusing that in one place, I needed to insert a space in *q/ *p to avoid accidentally starting a comment.)

Try it online!

Answer 3

Wolfram Language (Mathematica), 9 bytes

Resultant

Try it online!

Obligatory. And hat tip to OP...

Answer 4

J, 30 bytes

g~-/ .*@,g=.(+&-/}:@i.)&#{."{[

Try it online!

Straightforward but non-trivial answer.

Constructs the matrix and calculates the det -/ .*.

Answer 5

JavaScript (ES10), 171 bytes

Expects [p, q] where p and q are lists of coefficients, in descending order.

No built-in at all, so this is a bit lengthy.

a=>(D=m=>m+m?m.reduce((s,[v],i)=>s+v*(-1)**i*D(m.flatMap(([,...r])=>i--?[r]:[])),0):1)(a.flat().slice(2).map((_,y,v)=>v.map((_,x)=>~~a[i=1/a[1][y+1]?0:1][x-=i?y+~o:o=y])))

Try it online!

Building the Sylvester matrix

a.flat()                  // build a vector whose length is the sum of the
.slice(2)                 // lengths of p and q, minus 2
.map((_, y, v) =>         // for each entry at index y in this vector v[]:
  v.map((_, x) =>         //   for each entry at index x in v[]:
    ~~a[                  //     coerce to 0 if undefined:
      i =                 //       define i:
        1 / a[1][y + 1] ? //         if a[1][y + 1] is defined:
          0               //           set i = 0
        :                 //         else:
          1               //           set i = 1
    ][                    //
      x -=                //       subtract from x:
        i ?               //         if i is set:
          y + ~o          //           subtract y - o - 1
        :                 //         else:
          o = y           //           save y in o and subtract y
    ]                     //     read a[i][x]
  )                       //   end of inner map()
)                         // end of outer map()

Computing its determinant

The code is similar to the one used in my answer to this other challenge, except we make sure to return \$1\$ for the empty matrix.

Answer 6

Python3 + numpy, 155 bytes

from numpy import*
def f(p,q):P=len(p);Q=len(q);return linalg.det([[0]*i+p+[0]*(Q-2-i)for i in range(Q-1)]+[[0]*i+q+[0]*(P-2-i)for i in range(P-1)]or[[1]])

Try it online!

Answer 7

Excel, 164 bytes

Expects coefficient lists in descending order in vertical spilled ranges A1# and B1#.

=LET(
    a,A1#,
    b,B1#,
    c,ROWS(b),
    d,ROWS(a),
    f,LAMBDA(g,h,LET(i,SEQUENCE(,c+d-2)-SEQUENCE(h-1,,0),IFERROR(INDEX(g,IF(i,i,-1)),))),
    IFERROR(MDETERM(VSTACK(f(a,c),f(b,d))),1)
)

Answer 8

Haskell, 130 bytes

d=([]%);e%(f:g)=f!!0*d(tail<$>e++g)-(e++[f])%g;[]%_=1;_%_=0
x?[]=[x];x?(0:t)=[x++0:t]++(0:x)?t
u=drop 2.map(0*)
p#q=d$p?u q++q?u p

Try it online!

Determinant code stolen from xnor's https://codegolf.stackexchange.com/a/147820

Answer 9

MATL, 16 14 bytes

,&GnqXyY+]v&0|

Explanation

,       % Do twice
  &G    %   Push all inputs so far. The first time this implicitly takes p and
        %   pushes it. The second time this pushes p, then q
  nq    %   Number of elements minus 1
  Xy    %   Identity matrix of that size
  Y+    %   Two-dimensional convolution. The first time this implicitly takes q
        %   and uses it as first argument. The second time the first argument is p
]       % End
v       % Vertically concatenate the two matrices
&0|     % Determinant. Implicit display

Try it online! Or verify all test cases.

Answer 10

Vyxal, 27 25 bytes

@:∑2-:∇$-ʀ∇∆ZZƛ÷$¨VǓ;ÞfÞḊ

Try it Online!

-2 Thanks @The Thonnu!

A weird combination of Vyxal commands that unexpectedly gives the right result.
Welcome to golf it more!
Explanation:

@                           # Keep input and get lengths
 :∑2-:                      # Size of matrix 
       ∇$-ʀ∇               # Numbers of rotations
             ∆Z             # Pad with zeros to size
               Z            # Zip with numbers of rotations
                ƛ÷$¨VǓ;     # For both polynomials rotate every row
                       Þf   # Flatten
                         ÞḊ # Determinant

Answer 11

Charcoal, 81 bytes

≔⁺ＥΦηκ⁺Ｅκ⁰θＥΦθκ⁺Ｅκ⁰ηθ¿θ«≔⟦Ｅθ…⁺ιＥθ⁰Ｌθ⟧θＦθ¿⊖ＬιＦＥ⊟ιＥι×Φμ⁻πλ⎇ν¹×Ｘ±¹⁺Ｌιλκ⊞θκ⊞υ⊟⊟ιＩΣυ»1

Try it online! Link is to verbose version of code. Explanation:

≔⁺ＥΦηκ⁺Ｅκ⁰θＥΦθκ⁺Ｅκ⁰ηθ

Partly generate the Sylvester matrix.

¿θ«

If it's not empty then:

≔⟦Ｅθ…⁺ιＥθ⁰Ｌθ⟧θ

Pad the matrix so that it is square and wrap it in a list of matrix determinants to calculate.

Ｆθ¿⊖ＬιＦＥ⊟ιＥι×Φμ⁻πλ⎇ν¹×Ｘ±¹⁺Ｌιλκ⊞θκ⊞υ⊟⊟ιＩΣυ

Output the determinant as per my answer to Hankel transform of an integer sequence.

»1

Otherwise output the literal 1.

70 bytes by requiring the polynomials to be in ascending order of degree:

≔⊖ＥθＬιζＦ⌈ζＦθ⊞κ⁰≔…⁰Σζη⊞υ⟦⟧Ｆη≔ΣＥυＥ⁻⎇﹪λ²⮌ηηκ⁺κ⟦μ⟧υＩ↨Ｅ⮌υΠＥ駧θ÷λ⌈ζ⁻μ﹪λ⌈ζ±¹

Attempt This Online! Link is to verbose version of code. Takes input as a list of two polynomials in ascending order of degree. Explanation:

≔⊖ＥθＬιζ

Get the degrees of the polynomials.

Ｆ⌈ζＦθ⊞κ⁰

Pad them with zeros so that they can be readily cyclically indexed.

≔…⁰Σζη⊞υ⟦⟧Ｆη≔ΣＥυＥ⁻⎇﹪λ²⮌ηηκ⁺κ⟦μ⟧υ

Generate all of the permutations of [0..n+m) as per my answer to Hankel transform of an integer sequence.

Ｉ↨Ｅ⮌υΠＥ駧θ÷λ⌈ζ⁻μ﹪λ⌈ζ±¹

For each permutation, use the index for each element to determine the polynomial and offset and get the respective coefficient, then take the alternating sum of products, which is the resultant as required.

Answer 12

Jelly, 14 bytes

Ẉ’Ḷ0ẋ;€"ṚẎz0ÆḊ

Try it online!

Feels clumsy but it's at least shorter than I started with. Monadic link taking a list [q, p] (i.e. backwards).

   0ẋ             Repeat 0 to the length of each of
 ’Ḷ               the range from 0 to 2 less than
Ẉ                 the length of each coefficient list.
     ;            Append
       "Ṛ         the other coefficient list to each
      €           of the sublists.
         Ẏ        Concatenate the p rows and the q rows,
          z0      transpose with 0-padding (no effect on determinant),
            ÆḊ    and take the determinant.

Answer 13

Python + sympy, 30 bytes

from sympy import*
f=resultant

Attempt This Online!