Skip to content

Commit 1e98987

Browse files
committed
refactor: use RobustNonMonotoneLineSearch from LineSearch.jl
1 parent 5771fd0 commit 1e98987

File tree

4 files changed

+63
-180
lines changed

4 files changed

+63
-180
lines changed

src/NonlinearSolve.jl

Lines changed: 52 additions & 51 deletions
Original file line numberDiff line numberDiff line change
@@ -31,7 +31,7 @@ using LinearAlgebra: LinearAlgebra, ColumnNorm, Diagonal, I, LowerTriangular, Sy
3131
UpperTriangular, axpy!, cond, diag, diagind, dot, issuccess, istril,
3232
istriu, lu, mul!, norm, pinv, tril!, triu!
3333
using LineSearch: LineSearch, AbstractLineSearchAlgorithm, AbstractLineSearchCache,
34-
NoLineSearch
34+
NoLineSearch, RobustNonMonotoneLineSearch
3535
using LineSearches: LineSearches
3636
using LinearSolve: LinearSolve, LUFactorization, QRFactorization, ComposePreconditioner,
3737
InvPreconditioner, needs_concrete_A, AbstractFactorization,
@@ -105,54 +105,54 @@ include("algorithms/extension_algs.jl")
105105
include("utils.jl")
106106
include("default.jl")
107107

108-
# @setup_workload begin
109-
# nlfuncs = ((NonlinearFunction{false}((u, p) -> u .* u .- p), 0.1),
110-
# (NonlinearFunction{true}((du, u, p) -> du .= u .* u .- p), [0.1]))
111-
# probs_nls = NonlinearProblem[]
112-
# for (fn, u0) in nlfuncs
113-
# push!(probs_nls, NonlinearProblem(fn, u0, 2.0))
114-
# end
115-
116-
# nls_algs = (NewtonRaphson(), TrustRegion(), LevenbergMarquardt(),
117-
# PseudoTransient(), Broyden(), Klement(), DFSane(), nothing)
118-
119-
# probs_nlls = NonlinearLeastSquaresProblem[]
120-
# nlfuncs = ((NonlinearFunction{false}((u, p) -> (u .^ 2 .- p)[1:1]), [0.1, 0.0]),
121-
# (NonlinearFunction{false}((u, p) -> vcat(u .* u .- p, u .* u .- p)), [0.1, 0.1]),
122-
# (
123-
# NonlinearFunction{true}(
124-
# (du, u, p) -> du[1] = u[1] * u[1] - p, resid_prototype = zeros(1)),
125-
# [0.1, 0.0]),
126-
# (
127-
# NonlinearFunction{true}((du, u, p) -> du .= vcat(u .* u .- p, u .* u .- p),
128-
# resid_prototype = zeros(4)),
129-
# [0.1, 0.1]))
130-
# for (fn, u0) in nlfuncs
131-
# push!(probs_nlls, NonlinearLeastSquaresProblem(fn, u0, 2.0))
132-
# end
133-
134-
# nlls_algs = (LevenbergMarquardt(), GaussNewton(), TrustRegion(),
135-
# LevenbergMarquardt(; linsolve = LUFactorization()),
136-
# GaussNewton(; linsolve = LUFactorization()),
137-
# TrustRegion(; linsolve = LUFactorization()), nothing)
138-
139-
# @compile_workload begin
140-
# @sync begin
141-
# for T in (Float32, Float64), (fn, u0) in nlfuncs
142-
# Threads.@spawn NonlinearProblem(fn, T.(u0), T(2))
143-
# end
144-
# for (fn, u0) in nlfuncs
145-
# Threads.@spawn NonlinearLeastSquaresProblem(fn, u0, 2.0)
146-
# end
147-
# for prob in probs_nls, alg in nls_algs
148-
# Threads.@spawn solve(prob, alg; abstol = 1e-2, verbose = false)
149-
# end
150-
# for prob in probs_nlls, alg in nlls_algs
151-
# Threads.@spawn solve(prob, alg; abstol = 1e-2, verbose = false)
152-
# end
153-
# end
154-
# end
155-
# end
108+
@setup_workload begin
109+
nlfuncs = ((NonlinearFunction{false}((u, p) -> u .* u .- p), 0.1),
110+
(NonlinearFunction{true}((du, u, p) -> du .= u .* u .- p), [0.1]))
111+
probs_nls = NonlinearProblem[]
112+
for (fn, u0) in nlfuncs
113+
push!(probs_nls, NonlinearProblem(fn, u0, 2.0))
114+
end
115+
116+
nls_algs = (NewtonRaphson(), TrustRegion(), LevenbergMarquardt(),
117+
PseudoTransient(), Broyden(), Klement(), DFSane(), nothing)
118+
119+
probs_nlls = NonlinearLeastSquaresProblem[]
120+
nlfuncs = ((NonlinearFunction{false}((u, p) -> (u .^ 2 .- p)[1:1]), [0.1, 0.0]),
121+
(NonlinearFunction{false}((u, p) -> vcat(u .* u .- p, u .* u .- p)), [0.1, 0.1]),
122+
(
123+
NonlinearFunction{true}(
124+
(du, u, p) -> du[1] = u[1] * u[1] - p, resid_prototype = zeros(1)),
125+
[0.1, 0.0]),
126+
(
127+
NonlinearFunction{true}((du, u, p) -> du .= vcat(u .* u .- p, u .* u .- p),
128+
resid_prototype = zeros(4)),
129+
[0.1, 0.1]))
130+
for (fn, u0) in nlfuncs
131+
push!(probs_nlls, NonlinearLeastSquaresProblem(fn, u0, 2.0))
132+
end
133+
134+
nlls_algs = (LevenbergMarquardt(), GaussNewton(), TrustRegion(),
135+
LevenbergMarquardt(; linsolve = LUFactorization()),
136+
GaussNewton(; linsolve = LUFactorization()),
137+
TrustRegion(; linsolve = LUFactorization()), nothing)
138+
139+
@compile_workload begin
140+
@sync begin
141+
for T in (Float32, Float64), (fn, u0) in nlfuncs
142+
Threads.@spawn NonlinearProblem(fn, T.(u0), T(2))
143+
end
144+
for (fn, u0) in nlfuncs
145+
Threads.@spawn NonlinearLeastSquaresProblem(fn, u0, 2.0)
146+
end
147+
for prob in probs_nls, alg in nls_algs
148+
Threads.@spawn solve(prob, alg; abstol = 1e-2, verbose = false)
149+
end
150+
for prob in probs_nlls, alg in nlls_algs
151+
Threads.@spawn solve(prob, alg; abstol = 1e-2, verbose = false)
152+
end
153+
end
154+
end
155+
end
156156

157157
# Core Algorithms
158158
export NewtonRaphson, PseudoTransient, Klement, Broyden, LimitedMemoryBroyden, DFSane
@@ -172,8 +172,9 @@ export NewtonDescent, SteepestDescent, Dogleg, DampedNewtonDescent, GeodesicAcce
172172

173173
# Globalization
174174
## Line Search Algorithms
175-
export LineSearchesJL, NoLineSearch, RobustNonMonotoneLineSearch, LiFukushimaLineSearch
176-
export Static, HagerZhang, MoreThuente, StrongWolfe, BackTracking
175+
export LineSearchesJL, LiFukushimaLineSearch # FIXME: deprecated. use LineSearch.jl directly
176+
export Static, HagerZhang, MoreThuente, StrongWolfe, BackTracking # FIXME: deprecated
177+
export NoLineSearch, RobustNonMonotoneLineSearch
177178
## Trust Region Algorithms
178179
export RadiusUpdateSchemes
179180

src/algorithms/dfsane.jl

Lines changed: 3 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -19,9 +19,8 @@ For other keyword arguments, see [`RobustNonMonotoneLineSearch`](@ref).
1919
function DFSane(; σ_min = 1 // 10^10, σ_max = 1e10, σ_1 = 1, M::Int = 10, γ = 1 // 10^4,
2020
τ_min = 1 // 10, τ_max = 1 // 2, n_exp::Int = 2, max_inner_iterations::Int = 100,
2121
η_strategy::ETA = (fn_1, n, x_n, f_n) -> fn_1 / n^2) where {ETA}
22-
# linesearch = RobustNonMonotoneLineSearch(;
23-
# gamma = γ, sigma_1 = σ_1, M, tau_min = τ_min, tau_max = τ_max,
24-
# n_exp, η_strategy, maxiters = max_inner_iterations)
25-
linesearch = NoLineSearch()
22+
linesearch = RobustNonMonotoneLineSearch(;
23+
gamma = γ, sigma_1 = σ_1, M, tau_min = τ_min, tau_max = τ_max,
24+
n_exp, η_strategy, maxiters = max_inner_iterations)
2625
return GeneralizedDFSane{:DFSane}(linesearch, σ_min, σ_max, nothing)
2726
end

src/core/spectral_methods.jl

Lines changed: 5 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -118,7 +118,7 @@ end
118118
function SciMLBase.__init(prob::AbstractNonlinearProblem, alg::GeneralizedDFSane, args...;
119119
stats = empty_nlstats(), alias_u0 = false, maxiters = 1000,
120120
abstol = nothing, reltol = nothing, termination_condition = nothing,
121-
internalnorm::F = DEFAULT_NORM, maxtime = nothing, kwargs...) where {F}
121+
maxtime = nothing, kwargs...)
122122
timer = get_timer_output()
123123
@static_timeit timer "cache construction" begin
124124
u = __maybe_unaliased(prob.u0, alias_u0)
@@ -129,8 +129,7 @@ function SciMLBase.__init(prob::AbstractNonlinearProblem, alg::GeneralizedDFSane
129129
fu = evaluate_f(prob, u)
130130
@bb fu_cache = copy(fu)
131131

132-
linesearch_cache = __internal_init(prob, alg.linesearch, prob.f, fu, u, prob.p;
133-
stats, maxiters, internalnorm, kwargs...)
132+
linesearch_cache = init(prob, alg.linesearch, fu, u; stats, kwargs...)
134133

135134
abstol, reltol, tc_cache = init_termination_cache(
136135
prob, abstol, reltol, fu, u_cache, termination_condition)
@@ -166,7 +165,9 @@ function __step!(cache::GeneralizedDFSaneCache{iip};
166165
end
167166

168167
@static_timeit cache.timer "linesearch" begin
169-
linesearch_failed, α = __internal_solve!(cache.linesearch_cache, cache.u, cache.du)
168+
linesearch_sol = solve!(cache.linesearch_cache, cache.u, cache.du)
169+
linesearch_failed = !SciMLBase.successful_retcode(linesearch_sol.retcode)
170+
α = linesearch_sol.step_size
170171
end
171172

172173
if linesearch_failed

src/globalization/line_search.jl

Lines changed: 3 additions & 121 deletions
Original file line numberDiff line numberDiff line change
@@ -24,124 +24,6 @@ end
2424
Base.@deprecate LiFukushimaLineSearch(; nan_max_iter::Int = 5, kwargs...) LineSearch.LiFukushimaLineSearch(;
2525
nan_maxiters = nan_max_iter, kwargs...)
2626

27-
# """
28-
# RobustNonMonotoneLineSearch(; gamma = 1 // 10000, sigma_0 = 1, M::Int = 10,
29-
# tau_min = 1 // 10, tau_max = 1 // 2, n_exp::Int = 2, maxiters::Int = 100,
30-
# η_strategy = (fn₁, n, uₙ, fₙ) -> fn₁ / n^2)
31-
32-
# Robust NonMonotone Line Search is a derivative free line search method from DF Sane
33-
# [la2006spectral](@cite).
34-
35-
# ### Keyword Arguments
36-
37-
# - `M`: The monotonicity of the algorithm is determined by a this positive integer.
38-
# A value of 1 for `M` would result in strict monotonicity in the decrease of the L2-norm
39-
# of the function `f`. However, higher values allow for more flexibility in this reduction.
40-
# Despite this, the algorithm still ensures global convergence through the use of a
41-
# non-monotone line-search algorithm that adheres to the Grippo-Lampariello-Lucidi
42-
# condition. Values in the range of 5 to 20 are usually sufficient, but some cases may
43-
# call for a higher value of `M`. The default setting is 10.
44-
# - `gamma`: a parameter that influences if a proposed step will be accepted. Higher value
45-
# of `gamma` will make the algorithm more restrictive in accepting steps. Defaults to
46-
# `1e-4`.
47-
# - `tau_min`: if a step is rejected the new step size will get multiplied by factor, and
48-
# this parameter is the minimum value of that factor. Defaults to `0.1`.
49-
# - `tau_max`: if a step is rejected the new step size will get multiplied by factor, and
50-
# this parameter is the maximum value of that factor. Defaults to `0.5`.
51-
# - `n_exp`: the exponent of the loss, i.e. ``f_n=||F(x_n)||^{n\\_exp}``. The paper uses
52-
# `n_exp ∈ {1, 2}`. Defaults to `2`.
53-
# - `η_strategy`: function to determine the parameter `η`, which enables growth
54-
# of ``||f_n||^2``. Called as `η = η_strategy(fn_1, n, x_n, f_n)` with `fn_1` initialized
55-
# as ``fn_1=||f(x_1)||^{n\\_exp}``, `n` is the iteration number, `x_n` is the current
56-
# `x`-value and `f_n` the current residual. Should satisfy ``η > 0`` and ``∑ₖ ηₖ < ∞``.
57-
# Defaults to ``fn_1 / n^2``.
58-
# - `maxiters`: the maximum number of iterations allowed for the inner loop of the
59-
# algorithm. Defaults to `100`.
60-
# """
61-
# @kwdef @concrete struct RobustNonMonotoneLineSearch <:
62-
# AbstractNonlinearSolveLineSearchAlgorithm
63-
# gamma = 1 // 10000
64-
# sigma_1 = 1
65-
# M::Int = 10
66-
# tau_min = 1 // 10
67-
# tau_max = 1 // 2
68-
# n_exp::Int = 2
69-
# maxiters::Int = 100
70-
# η_strategy = (fn₁, n, uₙ, fₙ) -> fn₁ / n^2
71-
# end
72-
73-
# @concrete mutable struct RobustNonMonotoneLineSearchCache <:
74-
# AbstractNonlinearSolveLineSearchCache
75-
# f
76-
# p
77-
# ϕ
78-
# u_cache
79-
# fu_cache
80-
# internalnorm
81-
# maxiters::Int
82-
# history
83-
# γ
84-
# σ₁
85-
# M::Int
86-
# τ_min
87-
# τ_max
88-
# nsteps::Int
89-
# η_strategy
90-
# n_exp::Int
91-
# stats::NLStats
92-
# end
93-
94-
# function __internal_init(
95-
# prob::AbstractNonlinearProblem, alg::RobustNonMonotoneLineSearch, f::F, fu, u,
96-
# p, args...; stats, internalnorm::IN = DEFAULT_NORM, kwargs...) where {F, IN}
97-
# @bb u_cache = similar(u)
98-
# @bb fu_cache = similar(fu)
99-
# T = promote_type(eltype(fu), eltype(u))
100-
101-
# ϕ = @closure (f, p, u, du, α, u_cache, fu_cache) -> begin
102-
# @bb @. u_cache = u + α * du
103-
# fu_cache = evaluate_f!!(f, fu_cache, u_cache, p)
104-
# stats.nf += 1
105-
# return internalnorm(fu_cache)^alg.n_exp
106-
# end
107-
108-
# fn₁ = internalnorm(fu)^alg.n_exp
109-
# η_strategy = @closure (n, xₙ, fₙ) -> alg.η_strategy(fn₁, n, xₙ, fₙ)
110-
111-
# return RobustNonMonotoneLineSearchCache(
112-
# f, p, ϕ, u_cache, fu_cache, internalnorm, alg.maxiters,
113-
# fill(fn₁, alg.M), T(alg.gamma), T(alg.sigma_1), alg.M,
114-
# T(alg.tau_min), T(alg.tau_max), 0, η_strategy, alg.n_exp, stats)
115-
# end
116-
117-
# function __internal_solve!(cache::RobustNonMonotoneLineSearchCache, u, du; kwargs...)
118-
# T = promote_type(eltype(u), eltype(du))
119-
# ϕ = @closure α -> cache.ϕ(cache.f, cache.p, u, du, α, cache.u_cache, cache.fu_cache)
120-
# f_norm_old = ϕ(eltype(u)(0))
121-
# α₊, α₋ = T(cache.σ₁), T(cache.σ₁)
122-
# η = cache.η_strategy(cache.nsteps, u, f_norm_old)
123-
# f_bar = maximum(cache.history)
124-
125-
# for k in 1:(cache.maxiters)
126-
# f_norm = ϕ(α₊)
127-
# f_norm ≤ f_bar + η - cache.γ * α₊ * f_norm_old && return (false, α₊)
128-
129-
# α₊ *= clamp(α₊ * f_norm_old / (f_norm + (T(2) * α₊ - T(1)) * f_norm_old),
130-
# cache.τ_min, cache.τ_max)
131-
132-
# f_norm = ϕ(-α₋)
133-
# f_norm ≤ f_bar + η - cache.γ * α₋ * f_norm_old && return (false, -α₋)
134-
135-
# α₋ *= clamp(α₋ * f_norm_old / (f_norm + (T(2) * α₋ - T(1)) * f_norm_old),
136-
# cache.τ_min, cache.τ_max)
137-
# end
138-
139-
# return true, T(cache.σ₁)
140-
# end
141-
142-
# function callback_into_cache!(topcache, cache::RobustNonMonotoneLineSearchCache, args...)
143-
# fu = get_fu(topcache)
144-
# cache.history[mod1(cache.nsteps, cache.M)] = cache.internalnorm(fu)^cache.n_exp
145-
# cache.nsteps += 1
146-
# return
147-
# end
27+
function callback_into_cache!(topcache, cache::AbstractLineSearchCache, args...)
28+
LineSearch.callback_into_cache!(cache, get_fu(topcache))
29+
end

0 commit comments

Comments
 (0)