From 0056081d5acb2025795844166fd2218fe5c38a83 Mon Sep 17 00:00:00 2001 From: Benjamin Desef Date: Mon, 16 Mar 2026 17:02:30 +0100 Subject: [PATCH 1/9] New multithreaded API --- Project.toml | 4 +- src/MultistartOptimization.jl | 2 +- src/generic_api.jl | 8 ++- src/tiktak.jl | 105 ++++++++++++++++++++++++++-------- 4 files changed, 92 insertions(+), 27 deletions(-) diff --git a/Project.toml b/Project.toml index 71af462..817296a 100644 --- a/Project.toml +++ b/Project.toml @@ -1,11 +1,12 @@ name = "MultistartOptimization" uuid = "3933049c-43be-478e-a8bb-6e0f7fd53575" authors = ["Tamas K. Papp "] -version = "0.3.1" +version = "0.4" [deps] ArgCheck = "dce04be8-c92d-5529-be00-80e4d2c0e197" DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae" +OhMyThreads = "67456a42-1dca-4109-a031-0a68de7e3ad5" Sobol = "ed01d8cd-4d21-5b2a-85b4-cc3bdc58bad4" [weakdeps] @@ -18,5 +19,6 @@ MultistartOptimizationNLoptExt = "NLopt" ArgCheck = "1, 2" DocStringExtensions = "0.8, 0.9" NLopt = "0.5, 0.6, 1" +OhMyThreads = "0.8" Sobol = "1.3" julia = "1.10" diff --git a/src/MultistartOptimization.jl b/src/MultistartOptimization.jl index 2aa6169..142dbe6 100644 --- a/src/MultistartOptimization.jl +++ b/src/MultistartOptimization.jl @@ -3,9 +3,9 @@ module MultistartOptimization # NOTE: exports in included files using ArgCheck: @argcheck -using Base.Threads: @spawn, fetch using DocStringExtensions: FIELDS, FUNCTIONNAME, SIGNATURES, TYPEDEF using Sobol: SobolSeq, Sobol +using OhMyThreads include("generic_api.jl") include("tiktak.jl") diff --git a/src/generic_api.jl b/src/generic_api.jl index 3b90602..4176f06 100644 --- a/src/generic_api.jl +++ b/src/generic_api.jl @@ -58,9 +58,13 @@ function local_minimization(local_method, minimization_problem::MinimizationProb end """ -`$(FUNCTIONNAME)(multistart_method, local_method, minimization_problem)` +`$(FUNCTIONNAME)(multistart_method, local_method, minimization_problem, [scheduler]; + scheduler_options...)` -Multistart minimization using the given multistart and local methods. +Multistart minimization using the given multistart and local methods. The individual +minimizations are carried out using the `scheduler` from `OhMyThreads`. Instead of +passing a `Scheduler`, it is also possible to use `scheduler_options`, which will +automatically construct the most suitable `Scheduler`. """ function multistart_minimization end diff --git a/src/tiktak.jl b/src/tiktak.jl index 3f31a0a..a30ad7d 100644 --- a/src/tiktak.jl +++ b/src/tiktak.jl @@ -72,20 +72,15 @@ $(SIGNATURES) Evaluate and return points of an `N`-element Sobol sequence. -When `use_threads`, execution is parallelized using `Threads.@spawn`. +Execution is potentially parallelized using `tasks` Tasks. """ -function sobol_starting_points(minimization_problem::MinimizationProblem, N::Integer, - use_threads::Bool) +function sobol_starting_points(minimization_problem::MinimizationProblem, N::Integer, scheduler::Scheduler) (; objective, lower_bounds, upper_bounds) = minimization_problem s = SobolSeq(lower_bounds, upper_bounds) skip(s, N) # better uniformity - points = Iterators.take(s, N) + points = collect(Iterators.take(s, N)) _initial(x) = _objective_at_location(objective, x) - if use_threads - map(fetch, map(x -> @spawn(_initial(x)), points)) - else - map(_initial, points) - end + tmap(_initial, Base.promote_op(_initial, eltype(points)), points; scheduler) end """ @@ -93,41 +88,105 @@ $(SIGNATURES) Helper function to keep the `N` points with the lowest `value`. """ -function _keep_lowest(xs, N) +function _keep_lowest!(xs, N) @argcheck 1 ≤ N ≤ length(xs) - partialsort(xs, 1:N, by = p -> p.value) + partialsort!(xs, 1:N, by = p -> p.value) end +struct EnumeratedCatVector{T,X<:Tuple{AbstractVector{T},Vararg{AbstractVector{T}}}} <: AbstractVector{Tuple{Int,T}} + v::X +end + +Base.size(e::EnumeratedCatVector) = (sum(length, e.v, init=0),) +@generated function Base.getindex(e::EnumeratedCatVector, i::Int) + result = Expr(:block, Expr(:meta, :inline), Expr(:meta, :propagate_inbounds), + :(@boundscheck i ≤ 0 && throw(BoundsError(e, i))), + :(i_ = i) + ) + jmax = fieldcount(fieldtype(e, :v)) + for j in 1:fieldcount(fieldtype(e, :v))-1 + push!(result.args, :(len = length(e.v[$j])), :(if i ≤ len + return i_, @inbounds e.v[$j][begin+i-1] + else + i -= len + end)) + end + push!(result.args, :(return i_, e.v[$jmax][i])) + return result +end +Base.IndexStyle(::Type{X}) where {X<:EnumeratedCatVector} = IndexLinear() + """ $(SIGNATURES) Solve `minimization_problem` by using `local_method` within `multistart_method`. -When `use_threads`, initial point search is parallelized using `Threads.@spawn`. +Initial point search and subsequent minimizations are executed according to the scheduler +policy. The legacy argument `use_threads` enables default parallelization settings. More +control over multi-threading is possible by passing an explicit `scheduler` or, preferrably, +by using keyword arguments for constructing a scheduler from `OhMyThreads`. `prepend_points` should contain a vector of initial starting points that are prepended to the Sobol sequence. These are useful if a guess is available for the vicinity of the optimum. + +!!! warning + It is not advisable to use the `chunksize` option when constructing the scheduler, as the + construction of all possible initial point candidates and the actual minimization both use + parallelization, but have different numbers of points. It is recommended to use the keyword + interface instead of a `Scheduler` object, as different types of schedulers are optimal for + the different stages: a `StaticScheduler` for initial point selection and a `GreedyScheduler` + for the local minimization. """ function multistart_minimization(multistart_method::TikTak, local_method, - minimization_problem; - use_threads = true, - prepend_points = Vector{Vector{Float64}}()) + minimization_problem, scheduler::Union{Scheduler,Nothing}=nothing; + prepend_points = Vector{Vector{Float64}}(), + use_threads::Union{Nothing,Bool}=nothing, # legacy argument + ntasks::Union{Integer,Nothing}=nothing, + chunksize::Union{Integer,Nothing}=nothing, + minchunksize::Union{Integer,Nothing}=nothing, + chunking::Union{Bool,Nothing}=nothing, + nchunks::Union{Integer,Nothing}=nothing, + split::Union{OhMyThreads.Split,Symbol,Nothing}=nothing) + !isnothing(use_threads) && !isnothing(scheduler) && + error("Legacy argument use_threads cannot combined with new scheduler arguments") + !isnothing(scheduler) && !all(isnothing, (ntasks, chunksize, minchunksize, chunking, split)) && + error("Either a Scheduler or keyword options for the scheduler have to be specified") + if isnothing(use_threads) + use_threads = !(scheduler isa SerialScheduler || isnothing(ntasks) || isone(ntasks)) + end + if use_threads + if isnothing(scheduler) + initial_scheduler = StaticScheduler(; ntasks=isnothing(ntasks) ? Threads.nthreads() : ntasks, + chunksize, minchunksize, + chunking=isnothing(chunking) ? true : chunking, + split=isnothing(split) ? OhMyThreads.Consecutive() : split) + optimization_scheduler = GreedyScheduler(; ntasks=isnothing(ntasks) ? Threads.nthreads() : ntasks, + chunking=isnothing(chunking) ? false : chunking, + nchunks, chunksize, minchunksize, + split=isnothing(split) ? OhMyThreads.RoundRobin() : split) + else + initial_scheduler = optimization_scheduler = scheduler + end + else + initial_scheduler = optimization_scheduler = SerialScheduler() + end + # We now have at most two choices for the schedulers, so union splitting will work - this is type stable (; quasirandom_N, initial_N, θ_min, θ_max, θ_pow) = multistart_method (; objective, lower_bounds, upper_bounds) = minimization_problem @argcheck(all(x -> all(lower_bounds .≤ x .≤ upper_bounds), prepend_points), "prepend_points outside problem bounds") - quasirandom_points = sobol_starting_points(minimization_problem, quasirandom_N, - use_threads) - initial_points = _keep_lowest(quasirandom_points, initial_N) - all_points = vcat(map(x -> _objective_at_location(objective, x), prepend_points), - initial_points) - function _step(visited_minimum, (i, initial_point)) + quasirandom_points = sobol_starting_points(minimization_problem, quasirandom_N, initial_scheduler) + initial_points = _keep_lowest!(quasirandom_points, initial_N) + all_points = EnumeratedCatVector((map(Base.Fix1(_objective_at_location, objective), prepend_points), quasirandom_points)) + function _step((_, visited_minimum)::Tuple{Int,NamedTuple}, (i, initial_point)::Tuple{Int,NamedTuple}) θ = _weight_parameter(multistart_method, i) x = @. (1 - θ) * initial_point.location + θ * visited_minimum.location local_minimum = local_minimization(local_method, minimization_problem, x) local_minimum ≡ nothing && return visited_minimum - local_minimum.value < visited_minimum.value ? local_minimum : visited_minimum + # we don't care about the index parameter, but for type stability of the function, we'll always give back one + (0, local_minimum.value < visited_minimum.value ? (; local_minimum.location, local_minimum.value) : visited_minimum) end - foldl(_step, enumerate(Iterators.drop(all_points, 1)); init = first(all_points)) + # do not drop the first point - the local optimizer won't run on init! + treduce(_step, all_points; scheduler=optimization_scheduler, init=first(all_points))[2] end From c9f40c5918e61a83d60c5c33acb01fa99c0926a2 Mon Sep 17 00:00:00 2001 From: Benjamin Desef Date: Mon, 16 Mar 2026 17:26:59 +0100 Subject: [PATCH 2/9] Fix leftover from old approach --- src/tiktak.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/tiktak.jl b/src/tiktak.jl index a30ad7d..9d43258 100644 --- a/src/tiktak.jl +++ b/src/tiktak.jl @@ -72,7 +72,7 @@ $(SIGNATURES) Evaluate and return points of an `N`-element Sobol sequence. -Execution is potentially parallelized using `tasks` Tasks. +Execution is potentially parallelized using the strategy specified by `scheduler`. """ function sobol_starting_points(minimization_problem::MinimizationProblem, N::Integer, scheduler::Scheduler) (; objective, lower_bounds, upper_bounds) = minimization_problem From 76c0713efd0b3eda9c3bfa48f899f5e3d75c4fa0 Mon Sep 17 00:00:00 2001 From: Benjamin Desef Date: Wed, 8 Apr 2026 13:45:08 +0200 Subject: [PATCH 3/9] Expose both schedulers and define default scheduler functions --- src/tiktak.jl | 93 ++++++++++++++++++++++++++++----------------------- 1 file changed, 51 insertions(+), 42 deletions(-) diff --git a/src/tiktak.jl b/src/tiktak.jl index 9d43258..9fecee9 100644 --- a/src/tiktak.jl +++ b/src/tiktak.jl @@ -4,6 +4,8 @@ export TikTak +VERSION ≥ v"1.11.0-DEV.469" && eval(Expr(:public, :default_initial_point_scheduler, :default_local_scheduler)) + struct TikTak quasirandom_N::Int initial_N::Int @@ -93,6 +95,38 @@ function _keep_lowest!(xs, N) partialsort!(xs, 1:N, by = p -> p.value) end +""" +$(SIGNATURES) + +Constructs the scheduler used to compute the function values at the initial points. If +multiple threads are used, this is a `StaticScheduler`, to which arguments may be passed. +For a single-threaded run, this is a `SerialScheduler`. + +This function is not exported. +""" +default_initial_point_scheduler(; ntasks=Threads.nthreads(), chunksize=nothing, minchunksize=nothing, + split::Union{OhMyThreads.Split,Symbol}=OhMyThreads.Consecutive(), chunking::Bool=true) = + isone(Threads.nthreads()) || isone(ntasks) ? SerialScheduler() : + StaticScheduler(; ntasks, chunksize, minchunksize, split, chunking) + +""" +$(SIGNATURES) + +Constructs the scheduler used to invoke the local minimizer on the objective from various +initial points. If multiple threads are used, this is a `GreedyScheduler`, to which +arguments may be passed. For a single-threaded run, this is a `SerialScheduler`. + +This function is not exported. + +!!! warning + It is highly recommended to only use the `RoundRobin` splitting strategy in order to + resemble most closely the steps a serial execution would take. +""" +default_local_scheduler(; ntasks=Threads.nthreads(), nchunks=nothing, chunksize=nothing, minchunksize=nothing, + split::Union{OhMyThreads.Split,Symbol}=OhMyThreads.RoundRobin(), chunking::Bool=false) = + isone(Threads.nthreads()) || isone(ntasks) ? SerialScheduler() : + GreedyScheduler(; ntasks, nchunks, chunksize, minchunksize, split, chunking) + struct EnumeratedCatVector{T,X<:Tuple{AbstractVector{T},Vararg{AbstractVector{T}}}} <: AbstractVector{Tuple{Int,T}} v::X end @@ -122,61 +156,36 @@ $(SIGNATURES) Solve `minimization_problem` by using `local_method` within `multistart_method`. Initial point search and subsequent minimizations are executed according to the scheduler -policy. The legacy argument `use_threads` enables default parallelization settings. More -control over multi-threading is possible by passing an explicit `scheduler` or, preferrably, -by using keyword arguments for constructing a scheduler from `OhMyThreads`. +policy. The argument `use_threads` enables default parallelization settings; for more +fine-grained control, explicitly specify the schedulers for the two stages. `prepend_points` should contain a vector of initial starting points that are prepended to the Sobol sequence. These are useful if a guess is available for the vicinity of the optimum. !!! warning - It is not advisable to use the `chunksize` option when constructing the scheduler, as the - construction of all possible initial point candidates and the actual minimization both use - parallelization, but have different numbers of points. It is recommended to use the keyword - interface instead of a `Scheduler` object, as different types of schedulers are optimal for - the different stages: a `StaticScheduler` for initial point selection and a `GreedyScheduler` - for the local minimization. + If the local optimization method is deterministic, so is the multistart optimization, + provided the `local_scheduler` is serial (multi-threading initial point generation is + fine). For a multi-threaded `local_scheduler`, no guarantees can be made either to + visit the same starting points in each run. + While chunking can be enabled, it is strongly recommended to use the `RoundRobin` + method is to be preferred over the `Consecutive` splitting strategy for + `local_scheduler`, so that it is more likely to be mixing initial points from + probably better starting points. """ function multistart_minimization(multistart_method::TikTak, local_method, - minimization_problem, scheduler::Union{Scheduler,Nothing}=nothing; + minimization_problem; prepend_points = Vector{Vector{Float64}}(), - use_threads::Union{Nothing,Bool}=nothing, # legacy argument - ntasks::Union{Integer,Nothing}=nothing, - chunksize::Union{Integer,Nothing}=nothing, - minchunksize::Union{Integer,Nothing}=nothing, - chunking::Union{Bool,Nothing}=nothing, - nchunks::Union{Integer,Nothing}=nothing, - split::Union{OhMyThreads.Split,Symbol,Nothing}=nothing) - !isnothing(use_threads) && !isnothing(scheduler) && - error("Legacy argument use_threads cannot combined with new scheduler arguments") - !isnothing(scheduler) && !all(isnothing, (ntasks, chunksize, minchunksize, chunking, split)) && - error("Either a Scheduler or keyword options for the scheduler have to be specified") - if isnothing(use_threads) - use_threads = !(scheduler isa SerialScheduler || isnothing(ntasks) || isone(ntasks)) - end - if use_threads - if isnothing(scheduler) - initial_scheduler = StaticScheduler(; ntasks=isnothing(ntasks) ? Threads.nthreads() : ntasks, - chunksize, minchunksize, - chunking=isnothing(chunking) ? true : chunking, - split=isnothing(split) ? OhMyThreads.Consecutive() : split) - optimization_scheduler = GreedyScheduler(; ntasks=isnothing(ntasks) ? Threads.nthreads() : ntasks, - chunking=isnothing(chunking) ? false : chunking, - nchunks, chunksize, minchunksize, - split=isnothing(split) ? OhMyThreads.RoundRobin() : split) - else - initial_scheduler = optimization_scheduler = scheduler - end - else - initial_scheduler = optimization_scheduler = SerialScheduler() - end + use_threads::Bool=false, + initial_point_scheduler::Scheduler=use_threads ? default_initial_point_scheduler() : + SerialScheduler(), + local_scheduler::Scheduler=use_threads ? default_local_scheduler() : SerialScheduler()) # We now have at most two choices for the schedulers, so union splitting will work - this is type stable (; quasirandom_N, initial_N, θ_min, θ_max, θ_pow) = multistart_method (; objective, lower_bounds, upper_bounds) = minimization_problem @argcheck(all(x -> all(lower_bounds .≤ x .≤ upper_bounds), prepend_points), "prepend_points outside problem bounds") - quasirandom_points = sobol_starting_points(minimization_problem, quasirandom_N, initial_scheduler) + quasirandom_points = sobol_starting_points(minimization_problem, quasirandom_N, initial_point_scheduler) initial_points = _keep_lowest!(quasirandom_points, initial_N) all_points = EnumeratedCatVector((map(Base.Fix1(_objective_at_location, objective), prepend_points), quasirandom_points)) function _step((_, visited_minimum)::Tuple{Int,NamedTuple}, (i, initial_point)::Tuple{Int,NamedTuple}) @@ -188,5 +197,5 @@ function multistart_minimization(multistart_method::TikTak, local_method, (0, local_minimum.value < visited_minimum.value ? (; local_minimum.location, local_minimum.value) : visited_minimum) end # do not drop the first point - the local optimizer won't run on init! - treduce(_step, all_points; scheduler=optimization_scheduler, init=first(all_points))[2] + treduce(_step, all_points; scheduler=local_scheduler, init=first(all_points))[2] end From 05cd12bd72160f0d09ac75303e5744c658e22051 Mon Sep 17 00:00:00 2001 From: Benjamin Desef Date: Wed, 8 Apr 2026 13:45:22 +0200 Subject: [PATCH 4/9] Remove unnecessary allocation --- src/tiktak.jl | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/src/tiktak.jl b/src/tiktak.jl index 9fecee9..7206593 100644 --- a/src/tiktak.jl +++ b/src/tiktak.jl @@ -183,8 +183,10 @@ function multistart_minimization(multistart_method::TikTak, local_method, # We now have at most two choices for the schedulers, so union splitting will work - this is type stable (; quasirandom_N, initial_N, θ_min, θ_max, θ_pow) = multistart_method (; objective, lower_bounds, upper_bounds) = minimization_problem - @argcheck(all(x -> all(lower_bounds .≤ x .≤ upper_bounds), prepend_points), - "prepend_points outside problem bounds") + @argcheck(all(x -> begin + length(lower_bounds) == length(x) || throw(DimensionMismatch("prepend_points with invalid length")) + all(issorted, zip(lower_bounds, x, upper_bounds)) + end, prepend_points), "prepend_points outside problem bounds") quasirandom_points = sobol_starting_points(minimization_problem, quasirandom_N, initial_point_scheduler) initial_points = _keep_lowest!(quasirandom_points, initial_N) all_points = EnumeratedCatVector((map(Base.Fix1(_objective_at_location, objective), prepend_points), quasirandom_points)) From cc146354cf4624d50255f8ccd3d12c0e21656a73 Mon Sep 17 00:00:00 2001 From: Benjamin Desef Date: Wed, 8 Apr 2026 13:51:20 +0200 Subject: [PATCH 5/9] Split intransparent reduction function into two simpler ones: single-threaded folding and multi-threaded foreach --- src/tiktak.jl | 49 +++++++++++++++++++++++++++++++++++++++++++------ 1 file changed, 43 insertions(+), 6 deletions(-) diff --git a/src/tiktak.jl b/src/tiktak.jl index 7206593..887e4ae 100644 --- a/src/tiktak.jl +++ b/src/tiktak.jl @@ -189,15 +189,52 @@ function multistart_minimization(multistart_method::TikTak, local_method, end, prepend_points), "prepend_points outside problem bounds") quasirandom_points = sobol_starting_points(minimization_problem, quasirandom_N, initial_point_scheduler) initial_points = _keep_lowest!(quasirandom_points, initial_N) - all_points = EnumeratedCatVector((map(Base.Fix1(_objective_at_location, objective), prepend_points), quasirandom_points)) - function _step((_, visited_minimum)::Tuple{Int,NamedTuple}, (i, initial_point)::Tuple{Int,NamedTuple}) + _optimize( + multistart_method, + local_method, + minimization_problem, + EnumeratedCatVector((map(Base.Fix1(_objective_at_location, objective), prepend_points), initial_points)), + local_scheduler + ) +end + +function _optimize(multistart_method::TikTak, local_method, minimization_problem, all_points, ::SerialScheduler) + function _step(visited_minimum, (i, initial_point)) θ = _weight_parameter(multistart_method, i) x = @. (1 - θ) * initial_point.location + θ * visited_minimum.location local_minimum = local_minimization(local_method, minimization_problem, x) local_minimum ≡ nothing && return visited_minimum - # we don't care about the index parameter, but for type stability of the function, we'll always give back one - (0, local_minimum.value < visited_minimum.value ? (; local_minimum.location, local_minimum.value) : visited_minimum) + local_minimum.value < visited_minimum.value ? local_minimum : visited_minimum + end + foldl(_step, all_points; init=first(all_points)[2]) +end + +mutable struct VisitedMinimum{V,T} + location::V + value::T +end + +function _optimize(multistart_method::TikTak, local_method, minimization_problem, all_points, scheduler::Scheduler) + init = first(all_points)[2] + l = Threads.SpinLock() + visited_minimum = VisitedMinimum(init.location, init.value) + _step = let l=l, visited_minimum=visited_minimum + function _step((i, initial_point)::Tuple{Int,NamedTuple}) + θ = _weight_parameter(multistart_method, i) + x = @. (1 - θ) * initial_point.location + θ * visited_minimum.location + local_minimum = local_minimization(local_method, minimization_problem, x) + local_minimum ≡ nothing && return visited_minimum + if local_minimum.value < visited_minimum.value + lock(l) + if local_minimum.value < visited_minimum.value + visited_minimum.location = local_minimum.location + visited_minimum.value = local_minimum.value + end + unlock(l) + end + return + end end - # do not drop the first point - the local optimizer won't run on init! - treduce(_step, all_points; scheduler=local_scheduler, init=first(all_points))[2] + tforeach(_step, all_points; scheduler) + return (; location=visited_minimum.location, value=visited_minimum.value) end From f90b48a0b0b6869b57b7fabe579a15213adbcc3c Mon Sep 17 00:00:00 2001 From: Benjamin Desef Date: Wed, 8 Apr 2026 13:53:20 +0200 Subject: [PATCH 6/9] Avoid unnecessary allocations --- src/tiktak.jl | 38 +++++++++++++++++++++++++++++--------- 1 file changed, 29 insertions(+), 9 deletions(-) diff --git a/src/tiktak.jl b/src/tiktak.jl index 887e4ae..3afba1f 100644 --- a/src/tiktak.jl +++ b/src/tiktak.jl @@ -199,14 +199,32 @@ function multistart_minimization(multistart_method::TikTak, local_method, end function _optimize(multistart_method::TikTak, local_method, minimization_problem, all_points, ::SerialScheduler) - function _step(visited_minimum, (i, initial_point)) - θ = _weight_parameter(multistart_method, i) - x = @. (1 - θ) * initial_point.location + θ * visited_minimum.location - local_minimum = local_minimization(local_method, minimization_problem, x) - local_minimum ≡ nothing && return visited_minimum - local_minimum.value < visited_minimum.value ? local_minimum : visited_minimum + init = first(all_points)[2] + _step = let x=similar(init.location) + function _step(visited_minimum, (i, initial_point)) + θ = _weight_parameter(multistart_method, i) + @. x = (1 - θ) * initial_point.location + θ * visited_minimum.location + local_minimum = local_minimization(local_method, minimization_problem, x) + local_minimum ≡ nothing && return visited_minimum + local_minimum.value < visited_minimum.value ? local_minimum : visited_minimum + end + end + foldl(_step, all_points; init) +end + +@static if VERSION < v"1.12" + mutable struct OncePerTask{T, F} <: Function + const initializer::F + + OncePerTask{T}(initializer::Type{U}) where {T, U} = new{T,Type{U}}(initializer) + OncePerTask{T}(initializer::F) where {T, F} = new{T,F}(initializer) + OncePerTask{T,F}(initializer::F) where {T, F} = new{T,F}(initializer) + OncePerTask(initializer::Type{U}) where U = new{Base.promote_op(initializer), Type{U}}(initializer) + OncePerTask(initializer) = new{Base.promote_op(initializer), typeof(initializer)}(initializer) + end + @inline function (once::OncePerTask{T,F})() where {T,F} + get!(once.initializer, task_local_storage(), once)::T end - foldl(_step, all_points; init=first(all_points)[2]) end mutable struct VisitedMinimum{V,T} @@ -216,12 +234,14 @@ end function _optimize(multistart_method::TikTak, local_method, minimization_problem, all_points, scheduler::Scheduler) init = first(all_points)[2] + tls = OncePerTask{typeof(init.location)}(let x=init.location; () -> similar(x) end) l = Threads.SpinLock() visited_minimum = VisitedMinimum(init.location, init.value) - _step = let l=l, visited_minimum=visited_minimum + _step = let tls=tls, l=l, visited_minimum=visited_minimum function _step((i, initial_point)::Tuple{Int,NamedTuple}) θ = _weight_parameter(multistart_method, i) - x = @. (1 - θ) * initial_point.location + θ * visited_minimum.location + x = tls() + @. x = (1 - θ) * initial_point.location + θ * visited_minimum.location local_minimum = local_minimization(local_method, minimization_problem, x) local_minimum ≡ nothing && return visited_minimum if local_minimum.value < visited_minimum.value From edf9e4c466c1511e846fbcd78747144020a94e18 Mon Sep 17 00:00:00 2001 From: Benjamin Desef Date: Wed, 8 Apr 2026 13:57:35 +0200 Subject: [PATCH 7/9] Add comment why we define EnumeratedCatVector --- src/tiktak.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/tiktak.jl b/src/tiktak.jl index 3afba1f..81876de 100644 --- a/src/tiktak.jl +++ b/src/tiktak.jl @@ -127,6 +127,8 @@ default_local_scheduler(; ntasks=Threads.nthreads(), nchunks=nothing, chunksize= isone(Threads.nthreads()) || isone(ntasks) ? SerialScheduler() : GreedyScheduler(; ntasks, nchunks, chunksize, minchunksize, split, chunking) +# We need an indexable enumerate() and also want to iterate over the concatenation of two vectors without having to concatenate +# them. Define this simple struct to achieve both things. struct EnumeratedCatVector{T,X<:Tuple{AbstractVector{T},Vararg{AbstractVector{T}}}} <: AbstractVector{Tuple{Int,T}} v::X end From e6e45b31713b534bd15683a2136b26a6737f624d Mon Sep 17 00:00:00 2001 From: Benjamin Desef Date: Wed, 8 Apr 2026 14:16:00 +0200 Subject: [PATCH 8/9] Use TaskLocalValue from OMT instead of OncePerTask --- src/tiktak.jl | 21 +++------------------ 1 file changed, 3 insertions(+), 18 deletions(-) diff --git a/src/tiktak.jl b/src/tiktak.jl index 81876de..0d4073f 100644 --- a/src/tiktak.jl +++ b/src/tiktak.jl @@ -214,21 +214,6 @@ function _optimize(multistart_method::TikTak, local_method, minimization_problem foldl(_step, all_points; init) end -@static if VERSION < v"1.12" - mutable struct OncePerTask{T, F} <: Function - const initializer::F - - OncePerTask{T}(initializer::Type{U}) where {T, U} = new{T,Type{U}}(initializer) - OncePerTask{T}(initializer::F) where {T, F} = new{T,F}(initializer) - OncePerTask{T,F}(initializer::F) where {T, F} = new{T,F}(initializer) - OncePerTask(initializer::Type{U}) where U = new{Base.promote_op(initializer), Type{U}}(initializer) - OncePerTask(initializer) = new{Base.promote_op(initializer), typeof(initializer)}(initializer) - end - @inline function (once::OncePerTask{T,F})() where {T,F} - get!(once.initializer, task_local_storage(), once)::T - end -end - mutable struct VisitedMinimum{V,T} location::V value::T @@ -236,13 +221,13 @@ end function _optimize(multistart_method::TikTak, local_method, minimization_problem, all_points, scheduler::Scheduler) init = first(all_points)[2] - tls = OncePerTask{typeof(init.location)}(let x=init.location; () -> similar(x) end) + tlv = OhMyThreads.TaskLocalValue{typeof(init.location)}(let x=init.location; () -> similar(x) end) l = Threads.SpinLock() visited_minimum = VisitedMinimum(init.location, init.value) - _step = let tls=tls, l=l, visited_minimum=visited_minimum + _step = let tlv=tlv, l=l, visited_minimum=visited_minimum function _step((i, initial_point)::Tuple{Int,NamedTuple}) θ = _weight_parameter(multistart_method, i) - x = tls() + x = tlv[] @. x = (1 - θ) * initial_point.location + θ * visited_minimum.location local_minimum = local_minimization(local_method, minimization_problem, x) local_minimum ≡ nothing && return visited_minimum From b6c2e9a4da9374663e9147285f314b22110f7b11 Mon Sep 17 00:00:00 2001 From: Benjamin Desef Date: Thu, 9 Apr 2026 11:25:02 +0200 Subject: [PATCH 9/9] Again unify single- and multi-threaded approach to avoid code duplication Be careful with mutating optimizers, don't let them invalidate previous results --- src/tiktak.jl | 58 +++++++++++++++++++++++---------------------------- 1 file changed, 26 insertions(+), 32 deletions(-) diff --git a/src/tiktak.jl b/src/tiktak.jl index 0d4073f..c61fe0c 100644 --- a/src/tiktak.jl +++ b/src/tiktak.jl @@ -152,6 +152,17 @@ Base.size(e::EnumeratedCatVector) = (sum(length, e.v, init=0),) end Base.IndexStyle(::Type{X}) where {X<:EnumeratedCatVector} = IndexLinear() +struct NoLock <: Base.AbstractLock end +Base.lock(::NoLock) = nothing +Base.trylock(::NoLock) = true +Base.unlock(::NoLock) = nothing +Base.islocked(::NoLock) = false + +mutable struct VisitedMinimum{V,T} + location::V + value::T +end + """ $(SIGNATURES) @@ -191,49 +202,32 @@ function multistart_minimization(multistart_method::TikTak, local_method, end, prepend_points), "prepend_points outside problem bounds") quasirandom_points = sobol_starting_points(minimization_problem, quasirandom_N, initial_point_scheduler) initial_points = _keep_lowest!(quasirandom_points, initial_N) - _optimize( - multistart_method, - local_method, - minimization_problem, - EnumeratedCatVector((map(Base.Fix1(_objective_at_location, objective), prepend_points), initial_points)), - local_scheduler - ) -end + all_points = EnumeratedCatVector((map(Base.Fix1(_objective_at_location, objective), prepend_points), initial_points)) -function _optimize(multistart_method::TikTak, local_method, minimization_problem, all_points, ::SerialScheduler) init = first(all_points)[2] - _step = let x=similar(init.location) - function _step(visited_minimum, (i, initial_point)) - θ = _weight_parameter(multistart_method, i) - @. x = (1 - θ) * initial_point.location + θ * visited_minimum.location - local_minimum = local_minimization(local_method, minimization_problem, x) - local_minimum ≡ nothing && return visited_minimum - local_minimum.value < visited_minimum.value ? local_minimum : visited_minimum - end + if local_scheduler isa SerialScheduler + tlv = Ref(similar(init.location)) + l = NoLock() + else + tlv = OhMyThreads.TaskLocalValue{Base.RefValue{typeof(init.location)}}(let x=init.location; () -> Ref(similar(x)) end) + l = Threads.SpinLock() end - foldl(_step, all_points; init) -end - -mutable struct VisitedMinimum{V,T} - location::V - value::T -end - -function _optimize(multistart_method::TikTak, local_method, minimization_problem, all_points, scheduler::Scheduler) - init = first(all_points)[2] - tlv = OhMyThreads.TaskLocalValue{typeof(init.location)}(let x=init.location; () -> similar(x) end) - l = Threads.SpinLock() - visited_minimum = VisitedMinimum(init.location, init.value) + visited_minimum = VisitedMinimum(copy(init.location), init.value) _step = let tlv=tlv, l=l, visited_minimum=visited_minimum function _step((i, initial_point)::Tuple{Int,NamedTuple}) θ = _weight_parameter(multistart_method, i) - x = tlv[] + xref = local_scheduler isa SerialScheduler ? tlv : tlv[] + x = xref[] @. x = (1 - θ) * initial_point.location + θ * visited_minimum.location local_minimum = local_minimization(local_method, minimization_problem, x) local_minimum ≡ nothing && return visited_minimum if local_minimum.value < visited_minimum.value lock(l) if local_minimum.value < visited_minimum.value + # If local_minimization works in-place, local_minimum ≡ x - so the next iteration in this task would + # overwrite the minimum position. Hence, we swap our temporary with the previous minimum, which has now + # become obsolete. + xref[] = visited_minimum.location visited_minimum.location = local_minimum.location visited_minimum.value = local_minimum.value end @@ -242,6 +236,6 @@ function _optimize(multistart_method::TikTak, local_method, minimization_problem return end end - tforeach(_step, all_points; scheduler) + tforeach(_step, all_points; scheduler=local_scheduler) return (; location=visited_minimum.location, value=visited_minimum.value) end