The algorithm interface

This section starts a single, cohesive story that will weave through all documentation pages. We will incrementally build an iterative algorithm, enrich it with stopping criteria, and finally refine how it records (logs) its progress. Instead of presenting the API in the abstract, we anchor every concept in one concrete, tiny example you can copy & adapt.

Why an “interface” for algorithms? Iterative numerical methods nearly always share the same moving pieces:

immutable input (the mathematical problem you are solving),
immutable configuration (parameters and knobs of the chosen algorithm), and
mutable working data (current iterate, caches, diagnostics) that evolves as you step.

Bundling these together loosely without forcing one giant monolithic type makes it easier to:

reason about what is allowed to change and what must remain fixed,
write generic tooling (e.g. stopping logic, logging, benchmarking) that applies across many algorithms,
test algorithms in isolation by constructing minimal Problem/Algorithm pairs, and
extend behavior (add new stopping criteria, new logging events) without rewriting core loops.

The interface in this package formalizes those roles with three abstract types:

Problem: immutable, algorithm‑agnostic input data.
Algorithm: immutable configuration and parameters deciding how to iterate.
State: mutable data that evolves (current iterate, caches, counters, diagnostics).

It provides a framework for decomposing iterative methods into small, composable parts: concrete Problem/Algorithm/State types have to implement a minimal set of core functionality, and this package helps to stitch everything together and provide additional helper functionality such as stopping criteria and logging functionality.

Concrete example: Heron's method

To make everything tangible, we will work through a concrete example to illustrate the library's goals and concepts. Our running example is Heron's / Babylonian method for estimating $\sqrt{S}$. (see also the concise background on Wikipedia: Babylonian method (Heron's method)): Starting from an initial guess $x_0$, we may converge to the solution by iterating:

\[x_{k+1} = \frac{1}{2}\left(x_k + \frac{S}{x_k}\right)\]

We therefore suggest the following concrete implementations of the abstract types provided by this package: They are illustrative; various performance and generality questions will be left unaddressed to keep this example simple.

Algorithm types

using AlgorithmsInterface

struct SqrtProblem <: Problem
    S::Float64                # number whose square root we seek
end

struct HeronAlgorithm <: Algorithm
    stopping_criterion        # will be plugged in later (any StoppingCriterion)
end

mutable struct HeronState <: State
    iterate::Float64          # current iterate
    iteration::Int            # current iteration count
    stopping_criterion_state  # will be plugged in later (any StoppingCriterionState)
end

Initialization

In order to begin implementing the core parts of our algorithm, we start at the very beginning. There are two main entry points provided by the interface:

initialize_state constructs an entirely new state for the algorithm
initialize_state! (in-place) reset of an existing state.

An example implementation might look like:

function AlgorithmsInterface.initialize_state(
        problem::SqrtProblem, algorithm::HeronAlgorithm,
        stopping_criterion_state::StoppingCriterionState;
        kwargs...
    )
    x0 = rand()
    iteration = 0
    return HeronState(x0, iteration, stopping_criterion_state)
end

function AlgorithmsInterface.initialize_state!(
        problem::SqrtProblem, algorithm::HeronAlgorithm, state::HeronState;
        kwargs...
    )
    state.iteration = 0
    return state
end

Iteration steps

Algorithms define a mutable step via step!. For Heron's method:

function AlgorithmsInterface.step!(problem::SqrtProblem, algorithm::HeronAlgorithm, state::HeronState)
    S = problem.S
    x = state.iterate
    state.iterate = 0.5 * (x + S / x)
    return state
end

Note that we are only focussing on the actual algorithm, and not incrementing the iteration counter. These kinds of bookkeeping should be handled by the AlgorithmsInterface.increment! function, which will by default already increment the iteration counter. The following generic functionality is therefore enough for our purposes, and does not need to be defined. Nevertheless, if additional bookkeeping would be desired, this can be achieved by overloading that function:

function AlgorithmsInterface.increment!(state::State)
    state.iteration += 1
    return state
end

Running the algorithm

With these definitions in place you can already run (assuming you also choose a stopping criterion – added in the next section):

function heron_sqrt(x; maxiter = 10)
    prob = SqrtProblem(x)
    alg  = HeronAlgorithm(StopAfterIteration(maxiter))
    state = solve(prob, alg)  # allocates & runs
    return state.iterate
end

println("Approximate sqrt: ", heron_sqrt(16.0))

Approximate sqrt: 4.0

We will refine this example with better halting logic and logging shortly.

Reference: Core interface types & functions

Below are the automatic API docs for the core interface pieces. Read them after grasping the example above – the intent should now be clearer.

AlgorithmsInterface.initialize_state! — Method

state = initialize_state(
    problem::Problem, algorithm::Algorithm,
    [stopping_criterion_state::StoppingCriterionState];
    kwargs...
)
state = initialize_state!(
    problem::Problem, algorithm::Algorithm, state::State;
    kwargs...
)

Initialize a State based on a Problem and an Algorithm. The kwargs... should allow to initialize for example the initial point. This can be done in-place for state, then only values that did change have to be provided.

Note that since the returned state should also hold state.stopping_criterion_state, which will be used to keep the internal state of the stopping criterion, the out-of-place version receives this as one of its arguments. By default, that will be initialized separately through a call to initialize_stopping_state and provided as an argument.

On the other hand, the in-place version is not responsible for initializing the stopping_criterion_state, as that will be handled separately by initialize_stopping_state!.

Users that which to handle the stopping criterion initialization in initialize_state manually should overload the 2-argument version, while by default the 3-argument version should be implemented.