### Setup Roots.jl

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/reference.md

This setup code is used for running doctests and ensures the Roots package is available.

```julia
using Roots
```

--------------------------------

### Setup for Newton's Method Visualization

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/geometry-zero-finding.md

Imports necessary packages and defines a helper for computing the derivative using ForwardDiff. This setup is required for visualizing Newton's method.

```julia
using Roots
using Plots, ForwardDiff
Base.adjoint(f::Function)  = x  -> ForwardDiff.derivative(f, float(x)) # f' will compute derivative
```

--------------------------------

### Newton's Method Steps Visualization

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Visualizes the steps of Newton's method for the polynomial `x^5 - x - 1` starting from `0.1`. This example requires ForwardDiff for derivative calculation and plotting.

```julia
f(x) = x^5 - x - 1; # hide
D(f) = x -> ForwardDiff.derivative(f,float(x)) # hide
xs = [0.1] # x0
n = 15
for i in 1:(n-1) push!(xs, xs[end] - f(xs[end])/D(f)(xs[end])) end
ys = [zeros(Float64,n)\'map(f, xs)'][1:2n]
xs = xs[repeat(collect(1:n), inner=[2], outer=[1])]
plot(f, -1.25, 1.5, linewidth=3, legend=false)
plot!(zero, -1.25, 1.5, linewidth=3)
plot!(xs, ys)
#show(current())  # hide
savefig("newton.svg"); nothing # hide
```

--------------------------------

### Find Zero with Scalar and Secant Methods

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Illustrates using `find_zero` with a scalar initial guess, which defaults to `Order0()`, and explicitly using the secant method (`Order1()`). It also shows how to start the secant method with an iterable of two points.

```julia
julia> find_zero(f, 3) ≈ α₁  # find_zero(f, x0::Number) will use Order0()
true

julia> find_zero(f, 3, Order1()) ≈ α₁ # same answer, different method (secant)
true

julia> find_zero(f, (3, 2), Order1()) ≈ α₁ # start secant method with (3, f(3), (2, f(2))
true
```

--------------------------------

### Roots.jl Thukral8 Example

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Demonstrates a case where the Thukral8 algorithm in Roots.jl might produce a large, potentially misidentified zero due to relative tolerance issues with functions exhibiting sublinear growth.

```jldoctest
julia> find_zero(cbrt, 1, Roots.Thukral8())
1.725042287244107e23
```

--------------------------------

### Handle Convergence Failure with sin(x) at pi/2

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Demonstrates how to catch and handle convergence failures when using `find_zero` with `Order1()` for the `sin` function starting at `pi/2`.

```jldoctest
try  find_zero(sin, pi/2, Order1()) catch err  "Convergence failed" end
```

--------------------------------

### Solve root-finding problem with higher-order method

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Utilize higher-order methods by providing the function and its derivative to ZeroProblem. This example uses Roots.LithBoonkkampIJzerman to identify a zero at 0.0.

```julia
fx = ZeroProblem((f, D(f)), x0) # higher order methods can identify zero of this function

```

```julia
solve(fx, Roots.LithBoonkkampIJzerman(2,1), atol=0.0, rtol=0.0)

```

--------------------------------

### Initialize Chandrapatla State

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Implement the `init_state` method to set up the initial values for the Chandrapatla state object. This ensures the algorithm starts correctly, potentially with an initial bisection step.

```julia
function init_state(::Chandrapatla, F, x₀, x₁, fx₀, fx₁)
    a, b, fa, fb = x₁, x₀, fx₁, fx₀
    c, fc = a, fa
    ChandrapatlaState(b, a, c, fb, fa, fc)
end
```

--------------------------------

### SciPy-like Absolute Tolerance Example

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Illustrates a scenario where a conservative absolute tolerance, similar to SciPy's approach, fails to find a root for large values.

```jldoctest
julia> find_zero(x -> sqrt(eps()) - eps(x), (0,Inf))
9.981132799999999e7
```

--------------------------------

### Find zero using Bisection method

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Use `find_zero` with a bracket and `Bisection()` to find a root of a univariate function. The function `f(x) = cos(x) - x` is used as an example, and the output shows the found root and its corresponding function value.

```julia
julia> f(x) = cos(x) - x;

julia> x = find_zero(f, (0, pi/2), Bisection())
0.7390851332151607

julia> x, f(x)
(0.7390851332151607, 0.0)
```

--------------------------------

### Handle Convergence Failure with cbrt(x) at 1

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Illustrates convergence failure when using `find_zero` with `Order2()` for the cube root function `cbrt(x)` starting at `1`. This is due to a large second derivative.

```jldoctest
f(x) = cbrt(x);

x = try  find_zero(f, 1, Order2())  catch err  "Convergence failed" end	# all of 2, 5, 8, and 16 fail or diverge towards infinity
```

--------------------------------

### Secant Method and Bisection with Bracket

Source: https://context7.com/juliamath/roots.jl/llms.txt

Illustrates using the Secant method and Bisection method with bracket initial points for finding zeros.

```julia
Zsin = ZeroProblem(sin, (3, 4))
solve(Zsin, Secant())         # => 3.141592653589793
solve(Zsin, Order2())         # switch method in one word
solve(Zsin, Bisection(); xatol=1/16)   # with tolerance kwarg
```

--------------------------------

### Comparing `roots` from IntervalRootFinding with `find_zeros`

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

This snippet demonstrates the difference between rigorous interval-based root finding using `IntervalRootFinding.roots` and the approximate `find_zeros` from the Roots package. It shows how to extract numerical roots from the interval results.

```julia
julia> using IntervalArithmetic, IntervalRootFinding, Roots

julia> f(x) = sin(x) - 0.1*x^2 + 1
f (generic function with 1 method)

julia> rts = roots(f, -10..10)
4-element Vector{Root{Interval{Float64}}}:
 Root([3.14959, 3.1496], :unique)
 Root([-4.42654, -4.42653], :unique)
 Root([-3.10682, -3.10681], :unique)
 Root([-1.08205, -1.08204], :unique)

julia> find_zeros(f, -10, 10)
4-element Vector{Float64}:
 -4.426534982071949
 -3.1068165552293254
 -1.0820421327607177
  3.1495967624505226
```

```julia
julia> [find_zero(f, (interval(u).lo, interval(u).hi)) for u ∈ rts if u.status == :unique]
4-element Vector{Float64}:
  3.1495967624505226
 -4.426534982071949
 -3.1068165552293254
 -1.082042132760718
```

--------------------------------

### Find Zero with Bisection and Default Methods

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Demonstrates finding roots using Bisection and default methods. Bisection is used when a bracketing interval is specified or when the initial guess is not a scalar. The `A42` method is shown as an alternative for fewer function evaluations.

```julia
julia> using Roots

julia> f(x) = exp(x) - x^4;

julia> α₀, α₁, α₂ = -0.8155534188089607, 1.4296118247255556, 8.6131694564414;

julia> find_zero(f, (8,9), Bisection()) ≈ α₂ # a bisection method has the bracket specified
true

julia> find_zero(f, (-10, 0)) ≈ α₀ # Bisection is default if x in `find_zero(f, x)` is not scalar
true


julia> find_zero(f, (-10, 0), Roots.A42()) ≈ α₀ # fewer function evaluations than Bisection
true
```

--------------------------------

### Simplified Low-Overhead Root-Finding Functions

Source: https://context7.com/juliamath/roots.jl/llms.txt

Provides direct implementations of common root-finding algorithms (Secant, Bisection, Newton, DFree, Muller) for performance-critical scenarios, skipping convergence checks.

```julia
using Roots

f(x) = x^3 - 2x - 5

# Direct secant method
Roots.secant_method(f, 2.0, 2.1)          # => 2.0945514815423265

# Direct bisection
Roots.bisection(f, 1.0, 3.0)              # => 2.0945514815423265

# Direct Newton (pass f and f')
D(f) = x -> (f(x + 1e-8) - f(x)) / 1e-8
Roots.newton(f, D(f), 2.0)               # => 2.0945514815423265

# Derivative-free (Order0-like)
Roots.dfree(f, 2.0)

# Muller's method — uses three points, works for complex roots
Roots.muller(f, 1.0, 2.0, 3.0)
```

--------------------------------

### Find Zero with Default Method for cbrt(x)

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Shows that the default `find_zero` method successfully finds the root for `cbrt(x)` starting at `1` by identifying a bracket.

```jldoctest
x = find_zero(f, 1)
```

```jldoctest
x,  f(x)
```

--------------------------------

### Solving Zero Problems with Different Algorithms

Source: https://context7.com/juliamath/roots.jl/llms.txt

Demonstrates solving a `ZeroProblem` with default and specified algorithms. Non-convergence is handled by returning NaN.

```julia
Z = ZeroProblem(f, x0)
solve(Z)                              # default Order1 for scalar x0 => 0.7728829591492101
solve(Z, Order2(); atol=0.0, rtol=0.0)
```

```julia
f_bad(x) = cbrt(x) * exp(-x^2)
Zbad = ZeroProblem(f_bad, 0.1147)
solve(Zbad, Order1(), atol=0.0, rtol=0.0)   # => NaN
```

--------------------------------

### Simplified low-overhead functions

Source: https://context7.com/juliamath/roots.jl/llms.txt

Provides direct implementations of common root-finding algorithms like secant, bisection, and Newton, skipping convergence tracking for performance-critical applications.

```APIDOC
## Simplified low-overhead functions

### Description
For performance-critical inner loops the abstraction cost of `find_zero` may matter. These bare implementations skip convergence tracking and extra checks.

### Usage
```julia
Roots.algorithm(f, args...)
```

### Example
```julia
using Roots

f(x) = x^3 - 2x - 5

# Direct secant method
Roots.secant_method(f, 2.0, 2.1)          # => 2.0945514815423265

# Direct bisection
Roots.bisection(f, 1.0, 3.0)              # => 2.0945514815423265

# Direct Newton (pass f and f')
D(f) = x -> (f(x + 1e-8) - f(x)) / 1e-8
Roots.newton(f, D(f), 2.0)               # => 2.0945514815423265

# Derivative-free (Order0-like)
Roots.dfree(f, 2.0)

# Muller's method — uses three points, works for complex roots
Roots.muller(f, 1.0, 2.0, 3.0)
```
```

--------------------------------

### LithBoonkkampIJzerman Methods for Root Finding

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Demonstrates the `Roots.LithBoonkkampIJzerman` family of solvers, which generalize secant and Newton's methods by using a specified number of prior points (S) and derivatives (D). Higher values increase convergence rate at the cost of complexity.

```julia
find_zero(dfᵏs(f, 0), 2, Roots.LithBoonkkampIJzerman(3,0)) # like secant
```

```julia
find_zero(dfᵏs(f, 1), 2, Roots.LithBoonkkampIJzerman(2,1)) # like Newton
```

```julia
find_zero(dfᵏs(f, 2), 2, Roots.LithBoonkkampIJzerman(2,2)) # like Halley
```

--------------------------------

### MATLAB-Compatible `fzero` and `fzeros` Interface

Source: https://context7.com/juliamath/roots.jl/llms.txt

Presents a legacy interface (`fzero`, `fzeros`) for users familiar with MATLAB, supporting bracketing, scalar initial points, and finding all zeros within an interval.

```julia
using Roots

f(x) = exp(x) - x^4
α₀, α₁, α₂ = -0.8155534188089607, 1.4296118247255556, 8.6131694564414

fzero(f, 8, 9)    ≈ α₂   # bracketing
fzero(f, -10, 0)  ≈ α₀   # bracketing
fzero(f, 3)       ≈ α₁   # scalar initial point, Order0()
fzero(f, 3, order=1)      # Order1 (secant)
fzero(sin, big(3), order=16)  ≈ π

fzeros(f, -10, 10) ≈ [α₀, α₁, α₂]   # all zeros in interval
```

--------------------------------

### Find zero using A42 algorithm

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Illustrates using the `A42` algorithm directly with `find_zero` to find the root of `sin(x)` in the interval `(3, 4)`, which converges faster than the default bisection for `BigFloat`.

```julia
julia> find_zero(sin, (3,4), A42())
3.141592653589793
```

--------------------------------

### Define How Far Function for Projectile Range

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Defines 'howfar' to calculate the horizontal distance a projectile travels given a launch angle theta. It uses find_zero to solve for x where flight(x, theta) = 0, starting the search near the asymptote of the log term.

```jldoctest roots
julia> function howfar(theta)
         a = 200*cosd(theta)
         find_zero(x -> flight(x, theta), a-5)  # starting point has type determined by `theta`.
        end
howfar (generic function with 1 method)
```

--------------------------------

### Finding All Zeros of a Function

Source: https://context7.com/juliamath/roots.jl/llms.txt

Demonstrates using the `AllZeros` method to find all roots of a function within a given interval.

```julia
Zall = ZeroProblem(x -> exp(x) - x^4, (-10, 10))
solve(Zall, Roots.AllZeros())
```

--------------------------------

### Euler's Method Visualization

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/geometry-zero-finding.md

Visualizes Euler's method (irrational Halley) for root finding using a quadratic approximation. Requires `f`, `f'`, and `f''` to be defined.

```julia
L_f(x) = f(x) * f''(x) / (f'(x))^2

x0 = 1.4
x1 = x0 - 2 / (1 + sqrt(1 - 2L_f(x0))) * f(x0)/f'(x0)
t2(x) = f(x0) + f'(x0)*(x-x0) + f''(x0)/2 * (x - x0)^2

a, b = 1.1, 1.5
p = plot(f, a, b; legend=false, linewidth=3)
plot!(zero)
plot!(t2; color=:red, linewidth=3)

scatter!([x0, x1], [0,0]; markercolor=:blue)
annotate!([(x0,0,"x0", :bottom), (x1, 0, "x1", :bottom)])
scatter!([x0], [f(x0)]; markercolor=:blue)

scatter!([α],[0]; markercolor=:blue)
annotate!([(α, 0, "α", :top)])

p
```

--------------------------------

### Parameterized Functions and Broadcasting

Source: https://context7.com/juliamath/roots.jl/llms.txt

Shows how to solve `ZeroProblem` with parameterized functions `f(x, p)` and broadcasting over parameters.

```julia
g(x, p=2) = cos(x) - x/p
Zg = ZeroProblem(g, π/4)
solve(Zg, Order1())           # p=2 default => 1.0298665293222586
solve(Zg, Order1(), p=3)      # keyword     => 1.170120950002626
solve(Zg, Order1(), 4)        # positional  => 1.2523532340025887
```

```julia
ps = 2:5
[solve(Zg, Order1(), p) for p in ps]
```

--------------------------------

### Roots.AllZeros Documentation

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/reference.md

This snippet refers to the documentation for the AllZeros method in Roots.jl, used with `find_zeros` to heuristically scan for all zeros in an interval.

```julia
```@docs
Roots.AllZeros
```
```

--------------------------------

### Import Roots and ForwardDiff

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Import necessary packages for root finding and automatic differentiation.

```julia
using Roots, ForwardDiff
```

--------------------------------

### Define and Find Roots of a Function

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Define a function and use fzero to find its roots within specified intervals or with default settings. The bracketing method is shown, along with finding roots using a single point and higher-order methods.

```julia
julia> f(x) = exp(x) - x^4
f (generic function with 2 methods)
```

```julia
julia> fzero(f, 8, 9) ≈ α₂   # bracketing
true
```

```julia
julia> fzero(f, -10, 0) ≈ α₀
true
```

```julia
julia> fzeros(f, -10, 10) ≈ [α₀, α₁, α₂]
true
```

```julia
julia> fzero(f, 3) ≈ α₁      # default is Order0()
true
```

```julia
julia> fzero(sin, big(3), order=16)  ≈ π # uses higher order method
π
```

--------------------------------

### Find Zero Using Chandrapatla Solver

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Demonstrates how to use the custom Chandrapatla solver with the `find_zero` function for various mathematical functions.

```julia
find_zero(sin, (3,4), Chandrapatla())
find_zero(x -> exp(x) - x^4, (8,9), Chandrapatla())
find_zero(x -> x^5 - x - 1, (1,2), Chandrapatla())
```

--------------------------------

### find_zero with Order1 and zero tolerances

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Demonstrates find_zero with Order1 and zero tolerances, which results in a ConvergenceFailed error.

```julia
find_zero(f, x0, Roots.Order1(), atol=0.0, rtol=0.0) # error as no check on `|f(xn)|`

```

--------------------------------

### Chebyshev's Method Visualization

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/geometry-zero-finding.md

Visualizes Chebyshev's method for finding roots, showing the iterative steps and the inverse quadratic fit. Requires `x0`, `f`, `f'`, `L_f`, and `α` to be defined.

```julia
x1 = x0 - (1 + 1/2 * L_f(x0)) * f(x0) / f'(x0)
F(x, y) = -f''(x0)/(2f'(x0)^2) * (y-f(x0))^2 + y - f(x0)  - f'(x0) * (x- x0)

a, b = 1.1, 1.5
p = plot(f, a, b; legend=false, linewidth=3)
plot!(zero)
xs, ys = range(a, b, length=50), range(f(a), f(b), length=50);
zs = [F(x,y) for y ∈ ys, x ∈ xs];
contour!(xs, ys, zs; levels = [0], color=:red, linewidth=3)

scatter!([x0, x1], [0,0]; markermarkercolor=:blue)
annotate!([(x0,0,"x0", :bottom), (x1, 0, "x1", :bottom)])

scatter!([x0], [f(x0)]; markermarkercolor=:blue)

scatter!([α],[0]; markermarkercolor=:blue)
annotate!([(α, 0, "α", :top)])

p
```

--------------------------------

### Finding Critical Points via Derivative Zeros

Source: https://context7.com/juliamath/roots.jl/llms.txt

Demonstrates finding critical points of a function by finding the zeros of its derivative.

```julia
h(x) = 1/x^2 + x^3
find_zero(D(h), 1)    # => 0.9221079114817278  (critical point of h)
```

--------------------------------

### Find Zero with Order1 Method (Secant)

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Use the `Order1()` method for a more efficient algorithm than the default, suitable when speed is prioritized over robustness to initial guesses. Requires an initial guess.

```julia
julia> f(x) = 2x - exp(-x);

julia> x = find_zero(f, 1, Order1())
0.3517337112491958

julia> x, f(x)
(0.3517337112491958, -1.1102230246251565e-16)
```

--------------------------------

### Classical methods based on derivatives

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/reference.md

Documentation for classical root-finding methods that utilize function derivatives. These include Newton's and Halley's methods, among others.

```APIDOC
## Classical methods based on derivatives

We begin by describing the classical methods even though they are not necessarily recommended because they require more work of the user, as they give insight into why there are a variety of methods available.

The classical methods of [Newton](https://en.wikipedia.org/wiki/Newton%27s_method) and [Halley](https://en.wikipedia.org/wiki/Halley%27s_method) utilize information about the function and its derivative(s) in an iterative manner to converge to a zero of ``f(x)`` given an initial starting value.

Newton's method is easily described:

From an initial point, the next point in the iterative algorithm is found by identifying the intersection of the ``x`` axis with the tangent line of ``f`` at the initial point. This is repeated until convergence or the realization that convergence won't happen for the initial point. Mathematically,

``x_{n+1} = x_{n} - f(x_n)/f'(x_n).``

Some facts are helpful to understand the different methods available in `Roots`:

* For Newton's method there is a formula for the error: Set
  ``\epsilon_n = \alpha - x_n``, where ``\alpha`` is the zero, then
  ``\epsilon_{n+1} = -f''(\xi_n)/(2f'(\xi_n) \cdot \epsilon_n^2,``
  here ``\xi_n`` is some value between ``\alpha`` and ``x_n``.

* The error term, when of the form ``|\epsilon_{n+1}| \leq C\cdot|\epsilon_n|^2``, can be used to identify an interval around ``\alpha`` for which convergence is guaranteed. Such convergence is termed *quadratic* (order 2). For floating point solutions, quadratic convergence and a well chosen initial point can lead to convergence in 4 or 5 iterations. In general, convergence is termed order ``q`` when ``|\epsilon_{n+1}| \approx C\cdot|\epsilon_n|^q``

* The term ``-f''(\xi_n)/(2f'(\xi_n)`` indicates possible issues when ``f''`` is too big near ``\alpha`` or ``f'`` is too small near ``\alpha``. In particular if ``f'(\alpha) = 0``, there need not be quadratic convergence, and convergence can take many iterations. A zero for which ``f(x) = (x-\alpha)^{1+\beta}\cdot g(x)``, with ``g(\alpha) \neq 0`` is called *simple* when ``\beta=0`` and non-simple when ``\beta > 0``. Newton's method is quadratic near *simple zeros* and need not be quadratic near *non-simple* zeros.
As well, if ``f''`` is too big near ``\alpha``, or ``f'`` too small near ``\alpha``, or ``x_n`` too far from ``\alpha`` (that is, ``|\epsilon_n|>1``) the error might actually increase and convergence is not guaranteed.

* The explicit form of the error function can be used to guarantee convergence for functions with a certain shape (monotonic, convex functions where the sign of ``f''`` and ``f'`` don't change). Quadratic convergence may only occur once the algorithm is near the zero.

* The number of function evaluations per step for Newton's method is 2.

```@docs
Roots.Newton
Roots.Halley
Roots.QuadraticInverse
Roots.ChebyshevLike
Roots.SuperHalley
```

Newton and Halley's method are members of this family of methods:

```@docs
Roots.LithBoonkkampIJzerman{S,D}
```
```

--------------------------------

### Roots.Order0 Documentation

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/reference.md

This snippet refers to the documentation for the Order0 method in Roots.jl, which is a non-bracketing method.

```julia
```@docs
Roots.Order0
```
```

--------------------------------

### Numerically Find Inverse of x - sin(x)

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Demonstrates finding the inverse function of f(x) = x - sin(x) over the interval [0, 5*pi] and plotting it. Requires Plots and Roots packages.

```@example roots
using Plots, Roots;
f(x) = x - sin(x)
a, b = 0, 5pi
plot(inverse_function(f, a, b), f(a), f(b))
savefig("inversefunction.svg"); nothing # hide
```

--------------------------------

### Solve a root-finding problem with default settings

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Use the CommonSolve interface to define a root-finding problem and solve it using default settings. The default method is Order1.

```julia
f(x) = exp(-x) - x^3

```

```julia
x0 = 2.0

```

```julia
fx = ZeroProblem(f, x0)

```

```julia
solve(fx) ≈ 0.7728829591492101

```

--------------------------------

### fzero

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/reference.md

Provides a MATLAB-compatible interface for finding roots of a function. This interface is maintained for compatibility but is not recommended for new use.

```APIDOC
## fzero

### Description
Provides a MATLAB-compatible interface for finding roots of a function.

### Method
`fzero(f, x₀; options...)`

### Parameters
- `f`: The function for which to find the root.
- `x₀`: The initial guess or interval.
- `options...`: Additional keyword arguments for convergence criteria and tracking.

### Returns
- The approximate root of the function.
```

--------------------------------

### Solve a root-finding problem with Order2 method

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Solve a root-finding problem using the Order2 method with specified absolute and relative tolerances set to zero. This demonstrates passing specific methods and convergence criteria.

```julia
solve(fx, Order2(); atol=0.0, rtol=0.0) ≈ 0.7728829591492101

```

--------------------------------

### Find Zero with Unitful Units

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Demonstrates the integration of `Roots.jl` with the `Unitful` package, allowing root finding for functions with physical units.

```julia
julia> using Unitful

julia> s, m  = u"s", u"m";

julia> g, v₀, y₀ = 9.8*m/s^2, 10m/s, 16m;


julia> y(t) = -g*t^2 + v₀*t + y₀
y (generic function with 1 method)

julia> find_zero(y, 1s)  ≈ 1.886053370668014s
true
```

--------------------------------

### ZeroProblem and solve

Source: https://context7.com/juliamath/roots.jl/llms.txt

Define a root-finding problem once and solve it using different algorithms with various options. Handles non-convergence by returning NaN and supports parameterized functions.

```APIDOC
## ZeroProblem and solve

### Description
Define a root-finding problem and solve it using different algorithms. Supports parameterized functions and handles non-convergence.

### Usage
```julia
Z = ZeroProblem(f, x0)
solve(Z) # Uses default algorithm
solve(Z, algorithm; options...)
```

### Example
```julia
# Scalar function
f(x) = x^2 - 4
Z = ZeroProblem(f, 1.0)
solve(Z) # => 2.0
solve(Z, Order2(); atol=0.0, rtol=0.0)

# Parameterized function
g(x, p) = cos(x) - x/p
Zg = ZeroProblem(g, pi/4)
solve(Zg, Order1()) # p=2 default => 1.0298665293222586
solve(Zg, Order1(), p=3) # keyword argument => 1.170120950002626
solve(Zg, Order1(), 4) # positional argument => 1.2523532340025887

# Broadcasting over parameters
ps = 2:5
[solve(Zg, Order1(), p) for p in ps]

# Secant method with bracket initial points
Zsin = ZeroProblem(sin, (3, 4))
solve(Zsin, Secant()) # => 3.141592653589793
solve(Zsin, Order2()) # switch method
solve(Zsin, Bisection(); xatol=1/16) # with tolerance

# AllZeros method
Zall = ZeroProblem(x -> exp(x) - x^4, (-10, 10))
solve(Zall, Roots.AllZeros())
```
```

--------------------------------

### Solve with Bisection method and custom xatol

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Solve a root-finding problem using the Bisection method with a specified xatol (absolute tolerance for x). The Bisection method is guaranteed to converge.

```julia
fx = ZeroProblem(sin, (3,4))

```

```julia
solve(fx, Bisection(); xatol=1/16)

```

--------------------------------

### Find Zero with High-Order Convergence

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Demonstrates finding the root of `sin` using `BigFloat` and a high-order convergence method (`Order16()`), comparing its efficiency to `Order1()`.

```julia
julia> find_zero(sin, BigFloat(3.0), Order16()) ≈ π # 2 iterations to 6 using Order1()
true
```

--------------------------------

### Find Zero with Order2 Method (Steffensen)

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Employ the `Order2()` method for efficient root finding, particularly useful when a simple zero is not guaranteed. Requires an initial guess.

```julia
julia> f(x) = (x + 3) * (x - 1)^2;

julia> x = find_zero(f, -2, Order2())
-3.0

julia> x, f(x)
(-3.0, 0.0)
```

--------------------------------

### Problem-Algorithm-Solve Interface (Secant Method)

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Solves `sin(x) = 0` using the Secant method with the problem-algorithm-solve interface. A `ZeroProblem` is created with the function and initial interval, then solved with a `Secant` algorithm.

```julia
f(x) = sin(x);

x0 = (3, 4)

M = Secant()

Z = ZeroProblem(f, x0)

solve(Z, M)
```

--------------------------------

### `fzero` / `fzeros` — MATLAB-compatible interface

Source: https://context7.com/juliamath/roots.jl/llms.txt

A legacy interface for users familiar with MATLAB, providing `fzero` for finding a single root and `fzeros` for finding all roots within an interval.

```APIDOC
## `fzero` / `fzeros` — MATLAB-compatible interface

### Description
A legacy interface familiar to MATLAB users. Not recommended for new code, but maintained for compatibility.

### Usage
```julia
fzero(f, initial_guess, options...)
fzeros(f, interval)
```

### Example
```julia
using Roots

f(x) = exp(x) - x^4
α₀, α₁, α₂ = -0.8155534188089607, 1.4296118247255556, 8.6131694564414

fzero(f, 8, 9)    ≈ α₂   # bracketing
fzero(f, -10, 0)  ≈ α₀   # bracketing
fzero(f, 3)       ≈ α₁   # scalar initial point, Order0()
fzero(f, 3, order=1)      # Order1 (secant)
fzero(sin, big(3), order=16)  ≈ π

fzeros(f, -10, 10) ≈ [α₀, α₁, α₂]   # all zeros in interval
```
```

--------------------------------

### Root Finding with Automatic Differentiation (Halley's Method)

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Find roots using Halley's method with second derivatives computed automatically using ForwardDiff. The `D` and `dfᵏs` functions are used here.

```julia
find_zero((f, D(f), D(f,2)), 2, Roots.Halley())
```

--------------------------------

### Problem-Algorithm-Solve Interface (Order2 Method)

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Solves `sin(x) = 0` using the `Order2` method with the problem-algorithm-solve interface. This demonstrates changing the algorithm while keeping the same problem definition.

```julia
solve(Z, Order2())
```

--------------------------------

### Find zero of 1/x in (-1, 1)

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Demonstrates `find_zero` identifying a discontinuity at `0.0` for the function `1/x` within the bracket `(-1, 1)`, as the algorithm only looks for sign changes.

```julia
julia> find_zero(x -> 1/x, (-1, 1))
0.0
```

--------------------------------

### Root Finding with Automatic Differentiation (Newton's Method)

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Find roots using Newton's method with derivatives computed automatically using ForwardDiff. The derivative function `D` and helper `dfᵏs` are defined for this purpose.

```julia
function D(f, n::Int=1)
           n <= 0 && return f
           n == 1 && return x -> ForwardDiff.derivative(f,float(x))
           D(D(f,1),n-1)
       end

dfᵏs(f,k) = ntuple(i->D(f,i-1), Val(k+1)) # (f, f′, f′′, …)

find_zero((f,D(f)), 2, Roots.Newton())
```

--------------------------------

### Handling Non-Simple Zeros with Specialized Methods

Source: https://context7.com/juliamath/roots.jl/llms.txt

Compares standard methods with specialized ones like Order2B, Order1B, King, and Esser for functions with non-simple zeros (multiplicity > 1).

```julia
f(x) = (x + 3) * (x - 1)^2   # non-simple zero at x=1 (multiplicity 2)

# Standard Order2 needs ~51 function calls
find_zero(f, 2, Order2())       # converges but slowly

# Order2B: quadratic convergence independent of multiplicity — only 17 calls
find_zero(f, 2, Roots.Order2B())   # => ~1.0

# Order1B: superlinear, multiplicity-free
find_zero(f, 2, Roots.Order1B())

# King: superlinear derivative-free (phi ≈ 1.618), 2 function calls per step
find_zero(f, 2, Roots.King())

# Esser: order-2, derivative-free, 3 function calls per step
find_zero(f, 2, Roots.Esser())
```

--------------------------------

### Verify Secant Method Formula with SymPy

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/geometry-zero-finding.md

Verifies that the derived 'b' coefficient from the linear system matches the secant method formula by simplifying the difference to zero.

```julia
sm = x1 - y1 * (x1-x0)/(y1-y0)
simplify(sm - u[b])
```

--------------------------------

### Handle Convergence Failure for Polynomial Root

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Demonstrates convergence failure when using the default `find_zero` for the polynomial `x^5 - x - 1` with an initial guess of `0.1`, which is far from the actual root.

```jldoctest
f(x) = x^5 - x - 1;

x0 = 0.1

try find_zero(f, x0)  catch err  "Convergence failed" end
```

--------------------------------

### Parameterized function root finding

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Demonstrates finding roots of a parameterized function `g(x, p) = cos(x) - x/p` using `find_zero`. The parameter `p` can be passed as a keyword argument or a positional argument to solve different problems.

```julia
julia> g(x, p=1) = cos(x) - x/p;

julia> x0, M = (0, pi/2), Bisection()
((0, 1.5707963267948966), Bisection())

julia> find_zero(g, x0, M) # as before, solve cos(x) - x = 0 using default p=1
0.7390851332151607

julia> find_zero(g, x0, M; p=2) # solves cos(x) - x/2 = 0
1.0298665293222589

julia> find_zero(g, x0, M, 2) # positional argument; useful with broadcasting
1.0298665293222589
```

--------------------------------

### Find zero using bracketing method

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Use `find_zero` with a tuple representing an interval to find a root using a bracketing method. This method is guaranteed to converge for continuous functions.

```julia
find_zero(f, (-3/2,  -1/2)) ≈ -1.0983313019186336
true
```

--------------------------------

### Super-Halley Method Visualization

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/geometry-zero-finding.md

Visualizes the Super-Halley method, which is a specific case of a hyperbola-based root-finding approach. Requires `x0`, `f`, `f'`, `f''`, `L_f`, and `α` to be defined.

```julia
cn = -f'(x0)
bn = -f''(x0)/f'(x0)
an = -f''(x0)/(2f'(x0)^2) - bn/f'(x0)
x1 = x0 - (1 + 1/2 * (L_f(x0) / (1 + bn * f(x0) / f'(x0)))) * f(x0)/f'(x0)
F(x, y) = x - x0 + (y-f(x0)) * (1 + an * (y - f(x0))) / (bn * (y - f(x0)) + cn)


a, b= 1.1, 1.5
p = plot(f, a, b; legend=false, linewidth=3)
plot!(zero)
xs, ys = range(a, b, length=50), range(f(a), f(b), length=50);
zs = [F(x,y) for y ∈ ys, x ∈ xs];
contour!(xs, ys, zs; levels = [0], color=:red, linewidth=3)

scatter!([x0, x1], [0,0]; markercolor=:blue)
annotate!([(x0,0,"x0", :bottom), (x1, 0, "x1", :bottom)])
scatter!([x0], [f(x0)]; markercolor=:blue)

scatter!([α],[0]; markercolor=:blue)
annotate!([(α, 0, "α", :top)])

p
```

--------------------------------

### Use fzero with a real initial point

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

An alternative interface to Roots.jl, fzero(f, x0) calls a derivative-free method. The order parameter specifies the method order.

```julia
fzero(f, x0::Real; order=0)

```

--------------------------------

### Use fzero with a bracketing interval

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Use the fzero function with a bracketing interval [a, b]. This calls the find_zero algorithm with the Bisection method.

```julia
fzero(f, a::Real, b::Real)

```

--------------------------------

### Find Zero with Newton's Method and Automatic Differentiation

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Illustrates using Newton's method with automatic differentiation via the `ForwardDiff` package. A helper function `D` is defined to compute the derivative of a given function.

```julia
julia> using ForwardDiff

julia> D(f) = x -> ForwardDiff.derivative(f,float(x))
D (generic function with 1 method)
```

```julia
julia> f(x) = x^3 - 2x - 5
f (generic function with 1 method)

julia> x0 = 2
2

julia> find_zero((f, D(f)), x0, Roots.Newton()) ≈ 2.0945514815423265
true
```

--------------------------------

### find_zero with initial guess

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Use `find_zero` with a single initial guess to find a nearby zero using a non-bracketing method. These methods are often more efficient.

```APIDOC
## `find_zero` with initial guess

### Description
Finds a zero of a function near a given initial guess using a non-bracketing method.

### Method
`find_zero(f, initial_guess)`

### Parameters
#### Path Parameters
- **f** (Function) - The function for which to find a zero.
- **initial_guess** (Number) - An initial guess for the zero.

### Request Example
```julia
using Roots

f(x) =  x^5 - x + 1/2

find_zero(f,  0.6)
```

### Response
#### Success Response (Float64)
- Returns the approximate value of the zero found near the initial guess.

#### Response Example
```
0.550606579334135
```
```

--------------------------------

### fzeros

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/reference.md

Provides a MATLAB-compatible interface for finding roots of a function, potentially multiple roots. This interface is maintained for compatibility but is not recommended for new use.

```APIDOC
## fzeros

### Description
Provides a MATLAB-compatible interface for finding roots of a function, potentially multiple roots.

### Method
`fzeros(f, interval; options...)`

### Parameters
- `f`: The function for which to find the roots.
- `interval`: The interval to search for roots.
- `options...`: Additional keyword arguments for convergence criteria and tracking.

### Returns
- A vector of approximate roots of the function.
```

--------------------------------

### Define function with parameter

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Define a function `f(x, p)` where `p` is a scalar parameter. This is the first step in analyzing sensitivity.

```julia
f(x, p) = x^2 - p # p a scalar
p = 2
```

--------------------------------

### Find zeros for `sin(1/x)` near zero

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

For functions with infinitely many zeros near a point, like `sin(1/x)` near 0, `find_zeros` can identify a large number of them. Handle potential `NaN` results by defining the function appropriately.

```julia
julia> f(x) = iszero(x) ? NaN : sin(1/x);  # avoid sin(Inf) error

julia> rts = find_zeros(f, -1, 1);

julia> length(rts) # 88 zeros identified
88
```

--------------------------------

### Halley's Method for Root Finding

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Use Halley's method for cubic convergence. Requires the function, its first derivative, and its second derivative. These can be provided as a tuple.

```julia
fpp(x) = 6x;

x = Roots.find_zero((f, fp, fpp), 2, Roots.Halley())
x, f(x), sign(f(prevfloat(x)) * f(nextfloat(x)))
```

--------------------------------

### Find Zero using find_zero

Source: https://github.com/juliamath/roots.jl/blob/master/README.md

Demonstrates finding the zero of a function using `find_zero` with a specified interval. Ensure the `D` function and `med` are defined and the interval `(0, 1)` contains a root.

```julia
julia> find_zero(D(med), (0, 1)) ≈ median(as)
true
```

--------------------------------

### Chebyshev-like Methods

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/reference.md

Represents Chebyshev-like methods for finding zeros. These methods are part of a family that includes higher-order convergence properties.

```julia
using Roots

# Example usage:
f(x) = x^2 - 2
zero = find_zero(f, 1.0, ChebyshevLike())
println(zero)
```

--------------------------------

### Auto-differentiate with Bisection (pre-v1.9)

Source: https://github.com/juliamath/roots.jl/blob/master/docs/src/roots.md

Demonstrates an attempt to auto-differentiate with the Bisection method, which is not amenable to this approach. This highlights limitations in older Julia versions.

```julia
F(p) = find_zero(f, (zero(p), one(p)), Roots.Bisection(), p)
ForwardDiff.derivative(F, 1/2)
```