IntervalRootFinding.jl API
Externals
Types
Root
Object representing a possible root inside a given region. The field status
is either :unknown
or :unique
. If status
is :unique
then we know that there is a unique root of the function in question inside the given region.
Internally the status may also be :empty
for region guaranteed to contain no root, however such Root
s are discarded by default and thus never returned by the roots
function.
Fields
interval
: a region (eitherInterval
ofIntervalBox
) searched for roots.status
: the status of the region, valid values are:empty
,unkown
and:unique
.
Methods
roots(f, X, contractor=Newton, strategy=BreadthFirstSearch, tol=1e-15)
roots(f, deriv, X, contractor=Newton, strategy=BreadthFirstSearch, tol=1e-15)
roots(f, X, contractor, tol)
roots(f, deriv, X, contractor, tol)
Uses a generic branch and prune routine to find in principle all isolated roots of a function f:R^n → R^n
in a region X
, if the number of roots is finite.
Inputs:
f
: function whose roots will be foundX
:Interval
orIntervalBox
in which roots are searchedcontractor
: function that, when applied to the functionf
, determines the status of a given boxX
. It returns the new box and a symbol indicating the status. Current possible values areBisection
,Newton
andKrawczyk
deriv
: explicit derivative off
forNewton
andKrawczyk
strategy
:SearchStrategy
determining the order in which regions are processed.tol
: Absolute tolerance. If a region has a diameter smaller thantol
, it is returned with status:unknown
.
Function to solve a quadratic equation where the coefficients are intervals. Returns an array of intervals of the roots. Arguments a
, b
and c
are interval coefficients of x²
, x
and 1
respectively. The interval case differs from the non-interval case in that there might be three disjoint interval roots. In the third case, one interval root extends to −∞ and another extends to +∞. This algorithm finds the set of points where F.lo(x) ≥ 0
and the set of points where F.hi(x) ≤ 0
and takes the intersection of these two sets. Eldon Hansen and G. William Walster : Global Optimization Using Interval Analysis - Chapter 8
root_status(rt)
Return the status of a Root
.
Internals
Types
BBSearch{DATA}
Branch and bound search interface in element of type DATA.
This interface provide an iterable that perform the search.
There is currently three types of search supported BreadFirstBBSearch
, DepthFirstBBSearch
and KeyBBSearch
, each one processing the element of the tree in a different order. When subtyping one of these, the following methods must be implemented:
root_element(::BBSearch)
: return the element with which the search is startedprocess(::BBSearch, elem::DATA)
: return a symbol representing the action to perform with the elementelem
and an object of typeDATA
reprensenting the state of the element after processing (may returnelem
unchanged).bisect(::BBSearch, elem::DATA)
: return two elements of typeDATA
build by bisectingelem
Subtyping BBSearch
directly allows to have control over the order in which the elements are process. To do this the following methods must be implemented:
root_element(::BBSearch)
: return the first element to be processed. Use to build the initial tree.get_leaf_id!(::BBSearch, wt::BBTree)
: return the id of the next leaf that will be processed and remove it from the list of working leaves ofwt
.insert_leaf!(::BBSearch, wt::BBTree, leaf::BBLeaf)
: insert a leaf in the list of working leaves.
Valid symbols returned by the process function
:store
: the element is considered as final and is stored, it will not be further processed:bisect
: the element is bisected and each of the two resulting part will be processed:discard
: the element is discarded from the tree, allowing to free memory
BreadthFirstSearch <: BreadthFirstBBSearch
Type implementing the BreadthFirstBBSearch
interface for interval roots finding.
Fields
initial
: region (as aRoot
object) in which roots are searched.contractor
: contractor to use (Bisection
,Newton
orKrawczyk
)tol
: tolerance of the search
Contractor{F}
Abstract type for contractors.
DepthFirstSearch <: DepthFirstBBSearch
Type implementing the DepthFirstBBSearch
interface for interval roots finding.
Fields
initial
: region (as aRoot
object) in which roots are searched.contractor
: contractor to use (Bisection
,Newton
orKrawczyk
)tol
: tolerance of the search
Krawczyk{F, FP} <: Contractor{F}
Contractor type for the Krawczyk method.
Fields
f::F
: function whose roots are searchedf::FP
: derivative or jacobian off
(C::Krawczyk)(X, tol, α=where_bisect)
Contract an interval X
using Krawczyk operator and return the contracted interval together with its status.
Inputs
R
: Root object containing the interval to contract.tol
: Precision to which unique solutions are refined.α
: Point of bisection of intervals.
Newton{F, FP} <: Contractor{F}
Contractor type for the Newton method.
Fields
f::F
: function whose roots are searchedf::FP
: derivative or jacobian off
(C::Newton)(X, tol, α=where_bisect)
Contract an interval X
using Newton operator and return the contracted interval together with its status.
Inputs
R
: Root object containing the interval to contract.tol
: Precision to which unique solutions are refined.α
: Point of bisection of intervals.
Methods
branch_and_prune(X, contractor, strategy, tol)
Generic branch and prune routine for finding isolated roots using the given contractor to determine the status of a given box X
.
See the documentation of the roots
function for explanation of the other arguments.
data(leaf::BBLeaf)
Return the data stored in the leaf.
data(wt::BBTree)
Return all the data stored in a BBTree
as a list. The ordering of the elements is arbitrary.
Solves the system of linear equations using Gaussian Elimination. Preconditioning is used when the precondition
keyword argument is true
.
REF: Luc Jaulin et al., Applied Interval Analysis, pg. 72
Iteratively solves the system of interval linear equations and returns the solution set. Uses the Gauss-Seidel method (Hansen-Sengupta version) to solve the system. Keyword precondition
to turn preconditioning off. Eldon Hansen and G. William Walster : Global Optimization Using Interval Analysis - Chapter 5 - Page 115
newton1d
performs the interval Newton method on the given function f
with its derivative f′
and initial interval x
. Optional keyword arguments give the tolerances reltol
and abstol
. reltol
is the tolerance on the relative error whereas abstol
is the tolerance on |f(X)|, and a debug
boolean argument that prints out diagnostic information.
newton1d
performs the interval Newton method on the given function f
and initial interval x
. Optional keyword arguments give the tolerances reltol
and abstol
. reltol
is the tolerance on the relative error whereas abstol
is the tolerance on |f(X)|, and a debug
boolean argument that prints out diagnostic information.