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.00239210066773748+0.000685923025428412x^7-4.53097031258087\times {10}^{\mathrm{-6}}{x}^6-0.00164801014632734x^5+0.00111190210893585x^4+0.00232049175102826x^3-0.000814821101997568x^2-0.00159767495032094x,-0.000797538808884249-0.00120036529451165x^7-0.00177118364917771x^6+0.00150784758762038x^5+0.00376932816808507x^4+0.000171163088689716x^3-0.00139357959466411x^2+0.00145428191862992x,-0.964763821204430x+0.438529009807885+0.175411603893669x^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.249118969040x^8+473.172800945553x^7-2.81222933045108\times {10}^{\mathrm{-6}}{x}^6-0.0000147072040345022x^5+518.779744465642x^4-0.000370080042164914x^3-524.482546459755x^2+530.172317845435x-518.826092219543$

$\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.428067982309x^8-2.19177590934416\times {10}^{\mathrm{-7}}{x}^7-210.932454886683x^6+530.181566323192x^5+330.650841497779x^4-0.000159454935234903x^3+513.078143359541x^2-0.00399010296101806x+302.126536415563$

$\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.