Root finding

Updated 01 Nov 2026

•
Root finding
- •
  find $r$ such that $f(r)=0$ ; rewrite the equation in residual form first.
- •
  check continuity, domain restrictions, and possible multiple roots.
- •
  small $|f(x)|$ does not always imply small error in $x$ , especially when $f$ is flat. Near a simple root, $|x-r|\approx|f(x)|/|f'(r)|$ .
- •
  Method families
  - •
    Closed/bracketing -> Preserve an enclosing interval.
  - •
    Open -> Do not preserve a bracket.
  - •
    normally use one starting estimate for fixed point/Newton, two for secant.
•
Stopping tests
- •
  possible checks: bracket width/error bound, step change, residual size, iteration limit.
- •
  Tolerance scale
  - •
    Position -> atol, rtol.
  - •
    Residual -> ftol, in units of $f$ .
  - •
    compare the same stopping checks; equal numeric tolerances alone do not mean equal root accuracy.
- •
  Iteration record
  - •
    columns: $k,x_k,f(x_k)$ , step, and bracket width when present.
  - •
    count updates, not printed rows including $x_0$ .
  - •
    the code returns (root, updates, abs_residual).
  - •
    add logging inside the loop.
•
Bracketing methods
- •
  under continuity, a sign-changing interval contains at least one root.
  - •
    test endpoint roots first.
  - •
    avoid multiplying function values in software.
  - •
    $1/x$ changes sign across zero but is discontinuous and has no root.
  - •
    $x^2$ has a root without a sign change.
  - •
    If $f(c)=0$ , stop; otherwise replace the endpoint with the same sign as $f(c)$ , preserving opposite signs.
- •
  Bisection
  - •
    Candidate -> $c=(a+b)/2$ .
  - •
    halve the interval while preserving a sign change.
  - •
    reliable, predictable linear reduction of the bound; actual midpoint error need not halve exactly.
  - •
    After $n$ halvings
    
    •
    Width -> $w_n=(b_0-a_0)/2^n$ .
    
    •
    resulting-bracket midpoint error: at most $w_n/2$ .
  - •
    Candidate-count convention
    
    •
    $k=1$ uses the original bracket; evaluate each midpoint before updating its bracket.
    
    •
    Candidate error bound -> $(b_0-a_0)/2^k$ for candidate $k$ .
    
    •
    Required count -> $k\ge\lceil\log_2((b_0-a_0)/\varepsilon)\rceil$ , with at least one candidate, ensures midpoint error $\le\varepsilon$ absent an earlier exact-root exit.
  - •
    Example: $x^2-2$ , initial bracket $[1,2]$
    
    •
    candidate 1: $c_1=1.5$ , $f(c_1)=0.25$ ; retain $[1,1.5]$ .
    
    •
    candidate 2: $c_2=1.25$ , $f(c_2)=-0.4375$ ; retain $[1.25,1.5]$ .
    
    •
    candidate 3: $c_3=1.375$ , $f(c_3)=-0.109375$ .
- •
  Regula falsi
  - •
    Candidate -> $c=b-f(b)(b-a)/(f(b)-f(a))$ .
  - •
    use the line through endpoints while retaining the bracket.
  - •
    regula falsi is typically linear, not generically superlinear.
  - •
    one endpoint may remain fixed; the full bracket need not shrink to zero.
  - •
    Example: $x^2-2$ , initial bracket $[1,2]$
    
    •
    candidate 1: $c_1=4/3$ , $f(c_1)=-2/9$ .
    
    •
    candidate 2: $c_2=7/5$ , $f(c_2)=-1/25$ .
    
    •
    candidate 3: $c_3=24/17$ , $f(c_3)=-2/289$ .
    
    •
    right endpoint stays at 2.
  - •
    it does not halve the bracket. Use residual/step checks with the flat-function caution, or a certified bound if available.
- •
  Illinois
  - •
    halve a repeatedly retained endpoint's interpolation weight under the method's update rule to reduce regula-falsi stagnation.
  - •
    keep function values separately for signs and residuals.
- •
  Brent
  - •
    fast steps: Secant and inverse-quadratic interpolation with bracketing safeguards.
  - •
    fallback: Bisection when acceptance checks fail.
  - •
    Brent is not simply Newton plus bisection.
•
Fixed-point iteration
- •
  write the equation as $x=g(x)$ .
- •
  Iteration -> $x_{k+1}=g(x_k)$ .
- •
  Local derivative test
  - •
    Attractive -> $|g'(x^*)|<1$ is sufficient.
  - •
    Repelling -> $|g'(x^*)|>1$ .
  - •
    Inconclusive -> $|g'(x^*)|=1$ .
- •
  Interval guarantee
  - •
    $g$ maps a closed interval into itself and $|g'|\le q<1$ throughout it.
  - •
    contraction gives convergence and a unique fixed point there.
  - •
    Error bound -> $|x_k-x^*|\le q|x_k-x_{k-1}|/(1-q)$ .
- •
  convergence depends on rearrangement, root, and starting value.
- •
  quadratic example: $f(x)=x^2-3x+2$ :
  - •
    Rational rearrangement -> $x=-\frac{2}{x-3}$ , since $x(x-3)=-2$ .
    
    •
    Domain -> $x\ne3$ .
    
    •
    Derivative -> $g'=2/(x-3)^2$ .
    
    •
    root 1 is locally attractive.
    
    •
    root 2 is repelling.
  - •
    another rearrangement: $x=\frac{x^2+2}{3}$ .
    
    •
    at root 1, the derivative is $2/3$ , locally attractive.
    
    •
    at root 2, the derivative is $4/3$ , locally repelling.
  - •
    this second rearrangement is not universally divergent.
•
Newton–Raphson
- •
  Tangent update -> $x_{k+1}=x_k-\frac{f(x_k)}{f'(x_k)}$ .
- •
  Convergence
  - •
    at a simple root, convergence is rapid near the root when the derivative is usable.
  - •
    at a repeated root, ordinary Newton typically becomes linearly convergent.
- •
  possible failures: divergence, cycles, domain exit, or zero/very small $f'(x_k)$ .
  - •
    the derivative floor is scale-dependent, not universal.
- •
  quadratic example: $f(x)=x^2-3x+2$ :
  - •
    Derivative -> $f'(x)=2x-3$ .
  - •
    start $x_0=0$ : $x_1=2/3$ , $x_2=14/15$ , $x_3=254/255$ ; approaches root 1.
  - •
    start $x_0=1.5$ : Fails because $f'(1.5)=0$ .
•
Secant
- •
  uses two estimates to approximate the derivative, avoiding explicit evaluation of $f'(x)$ .
- •
  Update -> $x_{k+1}=x_k-\frac{f(x_k)(x_k-x_{k-1})}{f(x_k)-f(x_{k-1})}$ .
- •
  nonzero usable denominator; does not preserve a sign-changing bracket.
- •
  Evaluation cost
  - •
    Secant -> usually one new function evaluation per update after initialization.
  - •
    newton evaluates the function plus its derivative.
  - •
    fewest iterations is not necessarily cheapest.
- •
  Regula falsi also uses a secant line, but retains a bracket; open secant does not.
- •
•
Root-finding benchmark cases
- •
  Common cubic
  - •
    Function and derivative -> $f=x^3-x-2$ , $f'=3x^2-1$ .
  - •
    Bracket -> $[1,2]$ ; $f(1)=-2$ , $f(2)=4$ .
  - •
    Newton start -> $1.5$ .
  - •
    Secant starts -> $1,2$ .
  - •
    Fixed-point map/start -> $g(x)=\sqrt[3]{x+2}$ ; start at 1.
  - •
    Root -> About $1.521379706805$ .
  - •
    First updates
    
    •
    Bisection -> $1.5$ .
    
    •
    Regula falsi -> $4/3$ .
    
    •
    Newton -> $1.521739130435$ .
  - •
    for Newton, the module's residual-only $10^{-8}$ test needs 3 updates, not 4; the table includes $x_0$ .
  - •
    Module comparison at $10^{-6}$
    
    •
    Bisection -> 20 candidates for a position bound.
    
    •
    Regula falsi -> 13 updates for a residual bound.
    
    •
    unlike stopping criteria do not establish an equal-accuracy speed ranking.
- •
  Cosine phase
  - •
    Equation/map -> $\cos x=x$ , $g=\cos x$ .
  - •
    Root -> about $0.739085133215$ .
  - •
    on $[0,1]$ , $g$ maps into that interval and $|g'|\le\sin1<1$ .
  - •
    shows alternating linear convergence.
- •
  Transcendental practice
  - •
    Equation/bracket -> $e^x+x=0$ in $[-1,0]$ .
  - •
    comparison: bisection; newton start 0; secant starts -1 and 0.
  - •
    position bound: $10^{-10}$ needs 34 midpoint candidates unless an exact root occurs first.
- •
  Cycle practice
  - •
    function/start: Newton on $x^3-2x+2$ from 0.
  - •
    aequence: 0 -> 1 -> 0 -> 1.
  - •
    serivatives -2 and 1 are nonzero; a zero initial derivative is not the cause.
  - •
    report iteration-limit exhaustion.
- •
  Cumulative exercises
  - •
    solve taylor_sin(x, 15)-x/2=0 on $[1,\pi]$ .
  - •
    find the square-wave sum's first $0.5$ crossing.
  - •
    record the approximation order as well as root tolerance.

Python: five guarded solvers

•
guarded implementations of five classical root-finding algorithms with explicit stopping criteria.
•
Stopping behavior
- •
  endpoint/exact roots exit immediately
- •
  bisection: certified position bound
- •
  other four methods: step plus residual checks, unless the residual is exactly zero.
- •
  explicit reports for nonfinite values, unusable denominator/derivative, stagnation, and iteration-limit exhaustion.
•
Calls and defaults
- •
  demonstrations and assertions follow the definitions.
- •
  tolerances: atol=1e-10, rtol=0, ftol=1e-10; scale residual tolerance to the function.
- •
  iteration limit: 200.

import math

def root_limits(atol, rtol, ftol, max_iter):
      if not all(math.isfinite(v) for v in (atol, rtol, ftol)):
          raise ValueError("nonfinite tolerance")
      if atol <= 0 or rtol < 0 or ftol <= 0:
          raise ValueError("positive atol/ftol and nonnegative rtol required")
      if not isinstance(max_iter, int) or not 1 <= max_iter <= 100000:
          raise ValueError("invalid iteration limit")

def root_value(value):
      if not math.isfinite(value):
          raise ValueError("nonfinite iterate/function value")
      return value

def opposite(u, v):
      return (u < 0) != (v < 0)

def bracket(f, a, b):
      root_value(a); root_value(b)
      if a >= b:
          raise ValueError("need a < b")
      fa, fb = root_value(f(a)), root_value(f(b))
      if fa != 0 and fb != 0 and not opposite(fa, fb):
          raise ValueError("endpoints do not bracket a root")
      return fa, fb

def bisection(f, a, b, atol=1e-10, rtol=0.0, ftol=1e-10, max_iter=200):
      root_limits(atol, rtol, ftol, max_iter)
      fa, fb = bracket(f, a, b)
      if fa == 0: return a, 0, 0.0
      if fb == 0: return b, 0, 0.0
      for k in range(1, max_iter+1):
          c = root_value(a/2+b/2)
          fc = root_value(f(c))
          radius = b/2-a/2
          if fc == 0 or radius <= atol+rtol*abs(c):
              return c, k, abs(fc)
          if c == a or c == b:
              raise RuntimeError("bisection stagnated before position tolerance")
          if opposite(fa, fc): b, fb = c, fc
          else: a, fa = c, fc
      raise RuntimeError("bisection iteration limit")

def regula_falsi(f, a, b, atol=1e-10, rtol=0.0, ftol=1e-10, max_iter=200):
      root_limits(atol, rtol, ftol, max_iter)
      fa, fb = bracket(f, a, b)
      if fa == 0: return a, 0, 0.0
      if fb == 0: return b, 0, 0.0
      old = None
      for k in range(1, max_iter+1):
          scale = max(abs(fa), abs(fb))
          ua, ub = fa/scale, fb/scale
          denominator = ub-ua
          if denominator == 0:
              raise RuntimeError("regula falsi zero denominator")
          weight = -ua/denominator
          c = root_value((1-weight)*a+weight*b)
          fc = root_value(f(c))
          if fc == 0 or (old is not None and abs(c-old) <= atol+rtol*abs(c)
                         and abs(fc) <= ftol):
              return c, k, abs(fc)
          if c == a or c == b or c == old:
              raise RuntimeError("regula falsi stagnation")
          if opposite(fa, fc): b, fb = c, fc
          else: a, fa = c, fc
          old = c
      raise RuntimeError("regula falsi iteration limit")

def newton(f, df, x, atol=1e-10, rtol=0.0, ftol=1e-10,
             max_iter=200, derivative_floor=1e-14):
      root_limits(atol, rtol, ftol, max_iter)
      root_value(x)
      if not math.isfinite(derivative_floor) or derivative_floor <= 0:
          raise ValueError("positive finite derivative floor required")
      fx = root_value(f(x))
      if fx == 0: return x, 0, 0.0
      for k in range(1, max_iter+1):
          slope = root_value(df(x))
          if abs(slope) <= derivative_floor:
              raise RuntimeError("Newton derivative too small")
          y = root_value(x-fx/slope)
          fy = root_value(f(y))
          if fy == 0 or (abs(y-x) <= atol+rtol*abs(y) and abs(fy) <= ftol):
              return y, k, abs(fy)
          if y == x:
              raise RuntimeError("Newton stagnation")
          x, fx = y, fy
      raise RuntimeError("Newton iteration limit")

def fixed_point(g, f, x, atol=1e-10, rtol=0.0, ftol=1e-10, max_iter=200):
      root_limits(atol, rtol, ftol, max_iter)
      root_value(x)
      if root_value(f(x)) == 0: return x, 0, 0.0
      for k in range(1, max_iter+1):
          y = root_value(g(x))
          fy = root_value(f(y))
          if fy == 0 or (abs(y-x) <= atol+rtol*abs(y) and abs(fy) <= ftol):
              return y, k, abs(fy)
          if y == x:
              raise RuntimeError("fixed point stagnation with large residual")
          x = y
      raise RuntimeError("fixed point iteration limit")

def secant(f, x0, x1, atol=1e-10, rtol=0.0, ftol=1e-10, max_iter=200):
      root_limits(atol, rtol, ftol, max_iter)
      root_value(x0); root_value(x1)
      f0, f1 = root_value(f(x0)), root_value(f(x1))
      if f0 == 0: return x0, 0, 0.0
      if f1 == 0: return x1, 0, 0.0
      for k in range(1, max_iter+1):
          scale = max(abs(f0), abs(f1))
          u0, u1 = f0/scale, f1/scale
          denominator = u1-u0
          if abs(denominator) <= 1e-14:
              raise RuntimeError("secant denominator too small")
          y = root_value(x1-u1*(x1-x0)/denominator)
          fy = root_value(f(y))
          if fy == 0 or (abs(y-x1) <= atol+rtol*abs(y) and abs(fy) <= ftol):
              return y, k, abs(fy)
          if y == x1:
              raise RuntimeError("secant stagnation")
          x0, f0, x1, f1 = x1, f1, y, fy
      raise RuntimeError("secant iteration limit")

f = lambda x: x**3-x-2
df = lambda x: 3*x*x-1
def g(x):
      return math.copysign(abs(x+2)**(1/3), x+2)
answers = {
      "bisection": bisection(f, 1, 2),
      "regula falsi": regula_falsi(f, 1, 2),
      "Newton": newton(f, df, 1.5),
      "fixed point": fixed_point(g, f, 1),
      "secant": secant(f, 1, 2),
}
reference = 1.5213797068045676
for name, (root, updates, residual) in answers.items():
      assert abs(root-reference) <= 1.1e-10
      if name != "bisection": assert residual <= 1e-10
      print(f"{name}: root={root:.15f}, updates={updates}, residual={residual:.12g}")
assert bisection(lambda x: x, 0, 1) == (0, 0, 0.0)
assert newton(lambda x: x*x, lambda x: 2*x, 0) == (0, 0, 0.0)
root, n, residual = fixed_point(math.cos, lambda x: math.cos(x)-x, 1)
assert abs(root-0.7390851332151607) < 1e-10
print(f"cosine fixed point: root={root:.15f}, updates={n}, residual={residual:.12g}")
root, n, residual = bisection(lambda x: math.exp(x)+x, -1, 0)
assert n == 34 and abs(root+0.5671432904097838) < 1e-10
print(f"exp(x)+x: root={root:.15f}, candidates={n}, residual={residual:.12g}")
cycle_f = lambda x: x**3-2*x+2
cycle_df = lambda x: 3*x*x-2
z = 0.0
cycle = [z]
for _ in range(6):
      z -= cycle_f(z)/cycle_df(z)
      cycle.append(z)
assert cycle == [0, 1, 0, 1, 0, 1, 0]
failures = [
      (lambda: bisection(f, 2, 3), ValueError),
      (lambda: newton(lambda x: x*x-3*x+2, lambda x: 2*x-3, 1.5), RuntimeError),
      (lambda: secant(lambda x: x*x+1, -1, 1), RuntimeError),
      (lambda: newton(cycle_f, cycle_df, 0, max_iter=6), RuntimeError),
      (lambda: fixed_point(lambda x: x, f, 1), RuntimeError),
      (lambda: bisection(f, 1, 2, max_iter=1), RuntimeError),
      (lambda: bisection(lambda x: math.nan, 1, 2), ValueError),
      (lambda: bisection(f, 1, 2, atol=0), ValueError),
]
for call, expected in failures:
      try:
          call()
      except expected:
          pass
      else:
          raise AssertionError("expected root-finding failure was not reported")
print("root failure checks: 8 passed; Newton cycle:", cycle)

•

VBA

•
Guarded VBA implementations
- •
  method selector: RootSolve: "bisection", "regula", "newton", "fixed", or "secant".
- •
  shared validation; each actual algorithm is in the loop's Select Case.
•
Problem selector
- •
  "cubic" -> module 3's $x^3-x-2$ .
- •
  "quadratic" -> $x^2-3x+2$ .
- •
  "cosine", "exp", "cycle" -> the exercises above.
- •
  RootMap supports cubic/quadratic/cosine only; unsupported maps raise errors.
•
Calls and results
- •
  bisection: RootSolve("bisection","cubic",1,2).
- •
  newton: RootSolve("newton","cubic",1.5,0); the second starting estimate is unused for Newton/fixed point.
- •
  return value: variant array (root, update count, abs residual); assign it to result and access result(0).
- •
  bisection certifies position; others require step plus residual, except an exact zero.
- •
  optional arguments: atol, rtol, ftol, then maxIter.
- •
  derivative floor: $10^{-14}$ for these unit-scale examples.
•
Failure exercises
- •
  zero derivative: in DemoRoots, use result = RootSolve("newton", "quadratic", 1.5, 0); its handler reports the failure.
- •
  run DemoRootCycle for a caught iteration-limit error.
- •
  change the problem functions to reuse this template.

•

Private Function RootFunction(ByVal problem As String, ByVal x As Double) As Double
      Select Case problem
          Case "cubic": RootFunction = x ^ 3 - x - 2#
          Case "quadratic": RootFunction = x * x - 3# * x + 2#
          Case "cosine": RootFunction = Cos(x) - x
          Case "exp": RootFunction = Exp(x) + x
          Case "cycle": RootFunction = x ^ 3 - 2# * x + 2#
          Case Else: Err.Raise vbObjectError + 301, "RootFunction", "Unknown problem"
      End Select
End Function

Private Function RootDerivative(ByVal problem As String, ByVal x As Double) As Double
      Select Case problem
          Case "cubic": RootDerivative = 3# * x * x - 1#
          Case "quadratic": RootDerivative = 2# * x - 3#
          Case "cosine": RootDerivative = -Sin(x) - 1#
          Case "exp": RootDerivative = Exp(x) + 1#
          Case "cycle": RootDerivative = 3# * x * x - 2#
          Case Else: Err.Raise vbObjectError + 302, "RootDerivative", "Unknown problem"
      End Select
End Function

Private Function RootMap(ByVal problem As String, ByVal x As Double) As Double
      Select Case problem
          Case "cubic": RootMap = Sgn(x + 2#) * Abs(x + 2#) ^ (1# / 3#)
          Case "quadratic": RootMap = (x * x + 2#) / 3#
          Case "cosine": RootMap = Cos(x)
          Case Else: Err.Raise vbObjectError + 303, "RootMap", "No fixed-point map for this problem"
      End Select
End Function

Private Function RootOpposite(ByVal u As Double, ByVal v As Double) As Boolean
      RootOpposite = ((u < 0#) Xor (v < 0#))
End Function

Public Function RootSolve(ByVal method As String, ByVal problem As String, _
      ByVal first As Double, ByVal second As Double, _
      Optional ByVal atol As Double = 0.0000000001, _
      Optional ByVal rtol As Double = 0#, Optional ByVal ftol As Double = 0.0000000001, _
      Optional ByVal maxIter As Long = 200) As Variant
      Dim a As Double, b As Double, fa As Double, fb As Double
      Dim x As Double, fx As Double, c As Double, fc As Double, old As Double
      Dim scale As Double, ua As Double, ub As Double, den As Double, slope As Double
      Dim radius As Double, positionTol As Double, k As Long
      Dim closed As Boolean, haveOld As Boolean, converged As Boolean
      method = LCase$(method): problem = LCase$(problem)
      If atol <= 0# Or ftol <= 0# Or rtol < 0# Or maxIter < 1 Or maxIter > 100000 Then _
          Err.Raise vbObjectError + 304, "RootSolve", "Invalid tolerance or iteration limit"
      Select Case method
          Case "bisection", "regula": closed = True
          Case "newton", "fixed", "secant": closed = False
          Case Else: Err.Raise vbObjectError + 305, "RootSolve", "Unknown method"
      End Select
      If closed Or method = "secant" Then
          a = first: b = second
          If closed And a >= b Then Err.Raise vbObjectError + 306, "RootSolve", "Need a < b"
          fa = RootFunction(problem, a): fb = RootFunction(problem, b)
          If fa = 0# Then RootSolve = Array(a, 0, 0#): Exit Function
          If fb = 0# Then RootSolve = Array(b, 0, 0#): Exit Function
          If closed And Not RootOpposite(fa, fb) Then _
              Err.Raise vbObjectError + 307, "RootSolve", "Invalid bracket"
      Else
          x = first: fx = RootFunction(problem, x)
          If fx = 0# Then RootSolve = Array(x, 0, 0#): Exit Function
      End If
      For k = 1 To maxIter
          Select Case method
              Case "bisection"
                  c = a / 2# + b / 2#: radius = b / 2# - a / 2#
              Case "regula", "secant"
                  scale = Abs(fa): If Abs(fb) > scale Then scale = Abs(fb)
                  ua = fa / scale: ub = fb / scale: den = ub - ua
                  If den = 0# Then Err.Raise vbObjectError + 308, "RootSolve", "Zero line denominator"
                  If method = "secant" And Abs(den) <= 0.00000000000001 Then _
                      Err.Raise vbObjectError + 309, "RootSolve", "Secant denominator too small"
                  If method = "regula" Then
                      c = (1# + ua / den) * a - ua / den * b
                  Else
                      c = b - ub * (b - a) / den
                  End If
              Case "newton"
                  slope = RootDerivative(problem, x)
                  If Abs(slope) <= 0.00000000000001 Then _
                      Err.Raise vbObjectError + 310, "RootSolve", "Newton derivative too small"
                  c = x - fx / slope
              Case "fixed"
                  c = RootMap(problem, x)
          End Select
          fc = RootFunction(problem, c): positionTol = atol + rtol * Abs(c)
          converged = (fc = 0#)
          Select Case method
              Case "bisection"
                  If radius <= positionTol Then converged = True
              Case "regula"
                  If haveOld Then
                      If Abs(c - old) <= positionTol And Abs(fc) <= ftol Then converged = True
                  End If
              Case "secant"
                  If Abs(c - b) <= positionTol And Abs(fc) <= ftol Then converged = True
              Case Else
                  If Abs(c - x) <= positionTol And Abs(fc) <= ftol Then converged = True
          End Select
          Debug.Print method, k, c, fc
          If converged Then RootSolve = Array(c, k, Abs(fc)): Exit Function
          If closed Then
              If c = a Or c = b Then Err.Raise vbObjectError + 311, "RootSolve", "Bracket candidate stagnated"
              If method = "regula" And haveOld Then
                  If c = old Then Err.Raise vbObjectError + 312, "RootSolve", "Regula falsi stagnated"
              End If
              If RootOpposite(fa, fc) Then
                  b = c: fb = fc
              Else
                  a = c: fa = fc
              End If
          ElseIf method = "secant" Then
              If c = b Then Err.Raise vbObjectError + 313, "RootSolve", "Secant stagnated"
              a = b: fa = fb: b = c: fb = fc
          Else
              If c = x Then Err.Raise vbObjectError + 314, "RootSolve", "Open iteration stagnated"
              x = c: fx = fc
          End If
          old = c: haveOld = True
      Next k
      Err.Raise vbObjectError + 315, "RootSolve", "Iteration limit reached without convergence"
End Function

Public Sub DemoRoots()
      Dim method As Variant, result As Variant, first As Double, second As Double
      On Error GoTo Failed
      For Each method In Array("bisection", "regula", "newton", "fixed", "secant")
          first = 1#: second = 2#
          If method = "newton" Then first = 1.5
          result = RootSolve(CStr(method), "cubic", first, second)
          Debug.Print "RESULT", method, result(0), result(1), result(2)
      Next method
      result = RootSolve("fixed", "cosine", 1#, 0#)
      Debug.Print "cosine", result(0), result(1), result(2)
      Exit Sub
Failed:
      Debug.Print "FAILED: " & Err.Description
End Sub

Public Sub DemoRootCycle()
      Dim result As Variant
      On Error GoTo Failed
      result = RootSolve("newton", "cycle", 0#, 0#, 0.0000000001, 0#, 0.0000000001, 6)
      Debug.Print result(0)
      Exit Sub
Failed:
      Debug.Print "FAILED: " & Err.Description
End Sub

Linked references 1

numerical-methods-polban

◌ Explore connections in Graph view

Root finding

Root finding

Stopping tests

Bracketing methods

Fixed-point iteration

Newton–Raphson

Secant

Root-finding benchmark cases

Python: five guarded solvers

VBA

Linked references 1