Optimization.jl

The page is in the process of being translated.

There are some solvers that are available in the Optimization.jl package directly without the need to install any of the solver wrappers.

Methods

LBFGS: The popular quasi-Newton method that leverages limited memory BFGS approximation of the inverse of the Hessian. Through a wrapper over the L-BFGS-B fortran routine accessed from the LBFGSB.jl package. It directly supports box-constraints. This can also handle arbitrary non-linear constraints through a Augmented Lagrangian method with bounds constraints described in 17.4 of Numerical Optimization by Nocedal and Wright. Thus serving as a general-purpose nonlinear optimization solver available directly in Optimization.jl.

Examples

Unconstrained rosenbrock problem

using Optimization, OptimizationLBFGSB, Zygote

rosenbrock(x, p) = (p[1] - x[1])^2 + p[2] * (x[2] - x[1]^2)^2
x0 = zeros(2)
p = [1.0, 100.0]

optf = OptimizationFunction(rosenbrock, AutoZygote())
prob = Optimization.OptimizationProblem(optf, x0, p)
sol = solve(prob, LBFGSB())

retcode: Success
u: 2-element Vector{Float64}:
 0.9999997057368228
 0.999999398151528

With nonlinear and bounds constraints

function con2_c(res, x, p)
    res .= [x[1]^2 + x[2]^2, (x[2] * sin(x[1]) + x[1]) - 5]
end

optf = OptimizationFunction(rosenbrock, AutoZygote(), cons = con2_c)
prob = OptimizationProblem(optf, x0, p, lcons = [1.0, -Inf],
    ucons = [1.0, 0.0], lb = [-1.0, -1.0],
    ub = [1.0, 1.0])
res = solve(prob, LBFGSB(), maxiters = 100)

retcode: Success
u: 2-element Vector{Float64}:
 0.783397417853095
 0.6215211044097776

Train NN with Sophia

using OptimizationBase, OptimizationSophia, Lux, ADTypes, Zygote, MLUtils, Statistics, Random, ComponentArrays

x = rand(10000)
y = sin.(x)
data = MLUtils.DataLoader((x, y), batchsize = 100)

# Define the neural network
model = Chain(Dense(1, 32, tanh), Dense(32, 1))
ps, st = Lux.setup(Random.default_rng(), model)
ps_ca = ComponentArray(ps)
smodel = StatefulLuxLayer{true}(model, nothing, st)

function callback(state, l)
    state.iter % 25 == 1 && @show "Iteration: $(state.iter), Loss: $l"
    return l < 1e-1 ## Terminate if loss is small
end

function loss(ps, data)
    x_batch, y_batch = data
    ypred = [smodel([x_batch[i]], ps)[1] for i in eachindex(x_batch)]
    return sum(abs2, ypred .- y_batch)
end

optf = OptimizationFunction(loss, ADTypes.AutoZygote())
prob = OptimizationProblem(optf, ps_ca, data)

res = solve(prob, OptimizationSophia.Sophia(), callback = callback, epochs = 100)

retcode: Success
u: ComponentVector{Float32}(layer_1 = (weight = Float32[2.1733375; -1.1939192; … ; -1.9674702; -0.5265458;;], bias = Float32[0.034618918, -0.5436965, 0.24186288, 0.5109736, -0.42365417, -0.108434364, -0.6342865, -0.50556934, 0.31003064, -0.46368402  …  -0.68952954, 0.73083484, 0.0458446, 0.6992108, -0.27099517, 0.13010073, -0.019764228, -0.21060005, 0.4805976, -0.5327935]), layer_2 = (weight = Float32[-0.15994157 -0.26495138 … -0.30523074 0.014616169], bias = Float32[-0.073437415]))