Approximate a root of an odd degree polynomial

Question

Every odd degree polynomial has at least one real root. However this root does not have to be a rational number so your task is to output a sequence of rational numbers that approximates it.

Rules

• Your input is an odd degree polynomial \$f\$ with integer coefficients in any reasonable format.

• You must output a sequence of rational numbers \$(a_n)_{n=1}^{\infty}\$ such that \$f(\lim_{n\to\infty} a_n) =0\$ using the standard sequence rules.

• Your solution may fail because of finite precision however it must work in theory, assuming your language has arbitrary precision rational numbers. If your language doesn't support algebraic numbers, you may not assume that it does. In particular you may not assume that a find_roots built-in returns an infinite precision algebraic number if it doesn't.

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

Test cases

input -> example output

x^3-2 -> 1, 1.2, 1.25, 1.259, ...
x+7 -> -7, -7, -7, ...
x^5-x+1 -> -1.1, -1.167, -1.167304, ...
2x^7-3x^5-x^4-2x+1 -> 3, 0.7, 0.453 ...

Answer 1

05AB1E, 23 22 19 12 8 bytes

Lż{.Δβd

Try it online!

Takes as an input a number \$n\$ and the coefficients list, from the highest power to the smallest one, and returns \$a_n\$.

Requires the leading coefficient to be positive. If that is invalid, we need +2 bytes to fix it by adding *н after β

If this is invalid because we may not assume \$\tan\$ has infinite precision, then the following 13 bytes alternative could be used, which uses a spacing of \$\frac1 n\$ for the \$x\$s inside a bound:

*ÄO(D(Ÿ¹/.Δβd

Explanation

Finds the non-negative value of \$p(x)\$ (\$p(x) \cdot c_m\$ in the version which doesn't require a positive leading coefficient) with the smallest \$x\$ among the values \$ x = \tan(i) \$ for \$ 1 \leq i \leq n \$. The \$\tan\$ of integers is dense in \$\mathbb{R}\$ (Since \$ n\mod \pi\$ is dense in the range \$(0, \pi)\$ ), so we can get arbitrarily close to the target number, and since we always keep \$x\$s as \$n\$ increases the sequence is non-increasing. Additionally, uses the fact that base conversion is the same action as evaluating a polynomial at the base, and that 05AB1E's base conversion builtin works well with non-integers.

L          Get the range [1, 2, 3, ..., n]
ż         Calculate the tangent of each number in that range
{          And sort those
.Δ         Now find the first value such that
 β          When converting the polynomial from that base
 d          It is non-negative

Answer 2

Wolfram Language (Mathematica), 89 or 143 bytes

Answer #1

If it is enough to find one (any) real root we can use next facts:

• All real roots of odd degree polynomial are bounded inside interval \$[-M,M]\$ where \$M = \max\left(1,\sum_{k=0}^{2n}|a_k|\right)\$ (\$\sum_{k=0}^{2n}|a_k|\$ if coefficients are integers)
• Bisection method always converges to some root of continuous function if they exist in given interval

...and do not use extra-power built-ins.
Fortunately $MachineEpsilon works well with TIO. So next code produces stream of rational approximations for one of the real roots. Input: array of coefficients.

With[{c=#,m=Total[Abs/@#]},If[Abs[Echo@#1-#2]>$MachineEpsilon,p=Fold[x#1+#2&,0,c]/.x->#&;d=(#1+#2)/2;If[p@#1*p@d>0,#0[d,#2],#0[#1,d]]]&[-m,m]]&

Try it online!
Optimization is welcome!

Answer #2

The next code works well with Mathematica desktop (input is array of coefficients and output are approximations of all real roots):

Do[Echo[x/.NSolve[Fold[x#1+#2&,0,#],x,Reals,WorkingPrecision->n]],{n,∞}]&

But on TIO nor N, nor WorkingPrecision are not working(
So just like this (n>20 make no sense):

Do[Echo@Round[x/.NSolve[Fold[x#1+#2&,0,#],x,Reals],10.^-n],{n,20}]&

Try it online!

Answer 3

R, 91 90 bytes

\(x,n,a=polyroot(x),c=Re(a[!Im(a)])%/%1+0:n/n)c[order(sapply(c,\(y)x%*%y^seq(!x)/y)^2)][1]

Attempt This Online!

Uses the R polyroot function to find the (limited-precision) root a, then for each element n of the sequence, divides floor(a)-ceiling(a) into n steps and returns the one with the smallest absolute result after applying the polynomial.
The final result is still limited by R's precision*, but the approach should work in principle, and is not affected by the limited-precision of the polyroot function (as long as it is correct to the nearest integer).

The limit as n tends to infinity is the true real root of the polynomial; however, successive elements of the sequence may have increasing absolute differences from the true root, depending on whether partitioning the interval into n-1 or n steps allows a better solution. For a convergent-at-every-step sequence, divide into 2^n instead of n steps, for 4 more bytes: try it here.

---

(*) Note that this approach is actually very limited by R's precision, as the a[!Im(a)] part of the code (=get non-complex roots) fails quite easily due to floating-point inaccuracies: the test link overcomes this by redefining R's `!` (logical-NOT) operator to return TRUE as long as its argument is sufficiently small.

Answer 4

JavaScript (ES7), 74 bytes

-26 bytes (!) thanks to @CommandMaster

A port of @CommandMaster's method.

Expects (p)(n) where p is the coefficient list from lowest to highest.

p=>g=(n,m)=>n?g(n-1,p.reduce((q,c,i)=>q+c*(v=Math.tan(n))**i)<0|v>m?m:v):m

Try it online!

Commented

p =>                // p[] = list of coefficients
g = (               // g is a recursive function taking:
  n,                //   n = counter
  m                 //   m = minimum result, initially undefined
) =>                //
n ?                 // if n is not equal to 0:
  g(                //   do a recursive call:
    n - 1,          //     decrement n
    p.reduce(       //     using q as accumulator,
    (q, c, i) =>    //     for each coefficient c at index i:
      q +           //       add to q:
      c * (         //         c multiplied by v ** i
        v =         //         where v = tan(n)
        Math.tan(n) //
      ) ** i        //
    )               //     end of reduce()
    < 0 |           //     if the evaluation of the polynomial is negative
    v > m ?         //     or v is greater than m:
      m             //       keep m unchanged
    :               //     else:
      v             //       update m to v
  )                 //   end of recursive call
:                   // else:
  m                 //   stop and return m

Answer 5

Vyxal r, 8 bytes

∆P≬İ⌊=cḞ

Try it Online!

I think this is valid. Takes list of coefficients and then the nth term to generate.

Explained

∆P≬İ⌊=cḞ
∆P       # Get roots from polynomial with coefficients. Returns a list of exact algebraic numbers because sympy.
  ≬   c  # Get the first root where:
   İ⌊=   #   the root is real
       Ḟ # format to n decimal places

Answer 6

Charcoal, 30 bytes

⊞υΣ↔θ⊞υ±ΣυＦη«≔⊘Συζ§≔υ‹↨θζ⁰ζ⟦Ｉζ

Try it online! Takes the polynomial and a number of terms (but because it uses floats there's no point asking for more than 50). Assumes that the leading term is positive. Explanation: Simple binary search over the bound of roots given by @CommandMaster.

⊞υΣ↔θ

Set the upper bound to the sum of the absolute coefficients.

⊞υ±Συ

Set the lower bound to the negation of the upper bound.

Ｆη«

Repeat n times.

≔⊘Συζ

Calculate the average of the bounds.

§≔υ‹↨θζ⁰ζ

Evaluate the polynomial at that point and update the appropriate bound.

⟦Ｉζ

Also output the next approximation.