The QRGCD(f, g, x, eps) function returns univariate numeric polynomials u,v,d such that d is an approximate-GCD for the input polynomials (f, g). (See [1] for the definition of an approximate-GCD.)

With high probability, the output polynomials satisfy || u * f + v * g - d || < ||(f,g,u,v,d)|| * eps, || f - d * f1 || < ||f|| * eps, and || g - d * g1 || < ||g|| * eps, where the polynomials f1 and f2 are cofactors of f and g with respect to the divisor d.

The output polynomials u and v satisfy deg(u) < deg(g) - deg(d) and deg(v) < deg(f) - deg(d).

$\mathrm{with}\left(\mathrm{SNAP}\right)\&colon;$

$f≔\mathrm{expand}\left(\left(x-5\right)\left(x-\frac{1}{2}\right)\left(56x^8+83x^7+91x^4-92x^2+93x-91\right)\right)$

$f≔56x^{10}-225x^9+91x^6+\frac{271}{2}{x}^4+599x^3-\frac{1665}{2}{x}^2-\frac{633}{2}{x}^8-\frac{1001}{2}{x}^5+733x+\frac{415}{2}{x}^7-\frac{455}{2}$

$g≔\mathrm{expand}\left(\left(x-5\right)\left(x-\frac{1}{2}\right)\left(32x^8-37x^6+93x^5+58x^4+90x^2+53\right)\right)$

$g≔32x^{10}+43x^8+\frac{593}{2}{x}^7-546x^6+235x^4+278x^2-176x^9-\frac{173}{2}{x}^5-495x^3-\frac{583}{2}x+\frac{265}{2}$

$\mathrm{eps}≔1.0{10}^{-4}$

$\mathrm{eps}≔0.0001000000000$

$u,v,d≔\mathrm{QRGCD}\left(f\&comma;g\&comma;x\&comma;\mathrm{eps}\right)$

$u,v,d≔-0.00239210066773749+0.000685923025428401x^7-4.53097031259133\times {10}^{\mathrm{-6}}{x}^6-0.00164801014632735x^5+0.00111190210893583x^4+0.00232049175102823x^3-0.000814821101997611x^2-0.00159767495032097x,-0.000797538808884264-0.00120036529451162x^7-0.00177118364917765x^6+0.00150784758762042x^5+0.00376932816808507x^4+0.000171163088689713x^3-0.00139357959466413x^2+0.00145428191862992x,-0.964763821204428x+0.438529009807886+0.175411603893664x^2$

$\mathrm{evalb}\left(\mathrm{evalf}\left(\mathrm{norm}\left(\mathrm{expand}\left(uf+vg-d\right)\&comma;2\right)\right)<\mathrm{evalf}\left(\max \left(\mathrm{norm}\left(f\&comma;2\right)\&comma;\mathrm{norm}\left(g\&comma;2\right)\&comma;\mathrm{norm}\left(u\&comma;2\right)\&comma;\mathrm{norm}\left(v\&comma;2\right)\&comma;\mathrm{norm}\left(d\&comma;2\right)\right)\right)\mathrm{eps}\right)$

$\mathrm{true}$

$\mathrm{f1}≔\mathrm{Quotient}\left(f\&comma;d\&comma;x\right)$

$\mathrm{f1}≔319.249118969049x^8+473.172800945610x^7-2.81194595280338\times {10}^{\mathrm{-6}}{x}^6-0.0000147057886330958x^5+518.779744472732x^4-0.000370044593017768x^3-524.482546282526x^2+530.172318731588x-518.826087788793$

$\mathrm{evalb}\left(\mathrm{evalf}\left(\mathrm{norm}\left(\mathrm{expand}\left(f-\mathrm{f1}d\right)\&comma;2\right)\right)<\mathrm{evalf}\left(\mathrm{norm}\left(f\&comma;2\right)\right)\mathrm{eps}\right)$

$\mathrm{true}$

$\mathrm{g1}≔\mathrm{Quotient}\left(g\&comma;d\&comma;x\right)$

$\mathrm{g1}≔182.428067982314x^8-2.19151327159655\times {10}^{\mathrm{-7}}{x}^7-210.932454886558x^6+530.181566323828x^5+330.650841500969x^4-0.000159438983892106x^3+513.078143439310x^2-0.00398970411187071x+302.126538409815$

$\mathrm{evalb}\left(\mathrm{evalf}\left(\mathrm{norm}\left(\mathrm{expand}\left(g-\mathrm{g1}d\right)\&comma;2\right)\right)<\mathrm{evalf}\left(\mathrm{norm}\left(g\&comma;2\right)\right)\mathrm{eps}\right)$

$\mathrm{true}$

Corless, R. M.; Watt, S. M.; and Zhi, L. "QR Factoring to compute the GCD of Univariate Approximate Polynomials." IEEE Transactions on Signal Processing. Vol. 52(12), (2004): 3394-3402.