11"""
2- SimpleTrustRegion(; autodiff = AutoForwardDiff(), max_trust_radius::Real = 0.0,
3- initial_trust_radius::Real = 0.0, step_threshold::Real = 0.1,
4- shrink_threshold::Real = 0.25, expand_threshold::Real = 0.75,
5- shrink_factor::Real = 0.25, expand_factor::Real = 2.0, max_shrink_times::Int = 32)
2+ SimpleTrustRegion(; autodiff = AutoForwardDiff(), max_trust_radius = 0.0,
3+ initial_trust_radius = 0.0, step_threshold = nothing,
4+ shrink_threshold = nothing, expand_threshold = nothing,
5+ shrink_factor = 0.25, expand_factor = 2.0, max_shrink_times::Int = 32,
6+ nlsolve_update_rule = Val(false))
67
78A low-overhead implementation of a trust-region solver. This method is non-allocating on
89scalar and static array problems.
@@ -36,17 +37,22 @@ scalar and static array problems.
3637 `expand_threshold < r` (with `r` defined in `shrink_threshold`). Defaults to `2.0`.
3738 - `max_shrink_times`: the maximum number of times to shrink the trust region radius in a
3839 row, `max_shrink_times` is exceeded, the algorithm returns. Defaults to `32`.
40+ - `nlsolve_update_rule`: If set to `Val(true)`, updates the trust region radius using the
41+ update rule from NLSolve.jl. Defaults to `Val(false)`. If set to `Val(true)`, few of the
42+ radius update parameters -- `step_threshold = 0.05`, `expand_threshold = 0.9`, and
43+ `shrink_factor = 0.5` -- have different defaults.
3944"""
4045@kwdef @concrete struct SimpleTrustRegion <: AbstractNewtonAlgorithm
4146 autodiff = nothing
4247 max_trust_radius = 0.0
4348 initial_trust_radius = 0.0
4449 step_threshold = 0.0001
45- shrink_threshold = 0.25
46- expand_threshold = 0.75
47- shrink_factor = 0.25
50+ shrink_threshold = nothing
51+ expand_threshold = nothing
52+ shrink_factor = nothing
4853 expand_factor = 2.0
4954 max_shrink_times:: Int = 32
55+ nlsolve_update_rule = Val (false )
5056end
5157
5258function SciMLBase. __solve (prob:: NonlinearProblem , alg:: SimpleTrustRegion , args... ;
@@ -57,14 +63,27 @@ function SciMLBase.__solve(prob::NonlinearProblem, alg::SimpleTrustRegion, args.
5763 Δₘₐₓ = T (alg. max_trust_radius)
5864 Δ = T (alg. initial_trust_radius)
5965 η₁ = T (alg. step_threshold)
60- η₂ = T (alg. shrink_threshold)
61- η₃ = T (alg. expand_threshold)
62- t₁ = T (alg. shrink_factor)
66+ if alg. shrink_threshold === nothing
67+ η₂ = _unwrap_val (alg. nlsolve_update_rule) ? T (0.05 ) : T (0.25 )
68+ else
69+ η₂ = T (alg. shrink_threshold)
70+ end
71+ if alg. expand_threshold === nothing
72+ η₃ = _unwrap_val (alg. nlsolve_update_rule) ? T (0.9 ) : T (0.75 )
73+ else
74+ η₃ = T (alg. expand_threshold)
75+ end
76+ if alg. shrink_factor === nothing
77+ t₁ = _unwrap_val (alg. nlsolve_update_rule) ? T (0.5 ) : T (0.25 )
78+ else
79+ t₁ = T (alg. shrink_factor)
80+ end
6381 t₂ = T (alg. expand_factor)
6482 max_shrink_times = alg. max_shrink_times
6583 autodiff = __get_concrete_autodiff (prob, alg. autodiff)
6684
6785 fx = _get_fx (prob, x)
86+ norm_fx = norm (fx)
6887 @bb xo = copy (x)
6988 J, jac_cache = jacobian_cache (autodiff, prob. f, fx, x, prob. p)
7089 fx, ∇f = value_and_jacobian (autodiff, prob. f, fx, x, prob. p, jac_cache; J)
@@ -73,10 +92,17 @@ function SciMLBase.__solve(prob::NonlinearProblem, alg::SimpleTrustRegion, args.
7392 termination_condition)
7493
7594 # Set default trust region radius if not specified by user.
76- Δₘₐₓ == 0 && (Δₘₐₓ = max (norm (fx), maximum (x) - minimum (x)))
77- Δ == 0 && (Δ = Δₘₐₓ / 11 )
95+ Δₘₐₓ == 0 && (Δₘₐₓ = max (norm_fx, maximum (x) - minimum (x)))
96+ if Δ == 0
97+ if _unwrap_val (alg. nlsolve_update_rule)
98+ norm_x = norm (x)
99+ Δ = T (ifelse (norm_x > 0 , norm_x, 1 ))
100+ else
101+ Δ = T (Δₘₐₓ / 11 )
102+ end
103+ end
78104
79- fₖ = 0.5 * norm (fx) ^ 2
105+ fₖ = 0.5 * norm_fx ^ 2
80106 H = ∇f' * ∇f
81107 g = _restructure (x, ∇f' * _vec (fx))
82108 shrink_counter = 0
@@ -87,7 +113,7 @@ function SciMLBase.__solve(prob::NonlinearProblem, alg::SimpleTrustRegion, args.
87113 @bb Hδ = copy (x)
88114 dogleg_cache = (; δsd, δN_δsd, δN)
89115
90- for k in 1 : maxiters
116+ for _ in 1 : maxiters
91117 # Solve the trust region subproblem.
92118 δ = dogleg_method!! (dogleg_cache, ∇f, fx, g, Δ)
93119 @bb @. x = xo + δ
@@ -107,7 +133,7 @@ function SciMLBase.__solve(prob::NonlinearProblem, alg::SimpleTrustRegion, args.
107133 Δ = t₁ * Δ
108134 shrink_counter += 1
109135 shrink_counter > max_shrink_times && return build_solution (prob, alg, x, fx;
110- retcode = ReturnCode. ConvergenceFailure )
136+ retcode = ReturnCode. ShrinkThresholdExceeded )
111137 end
112138
113139 if r ≥ η₁
@@ -121,12 +147,22 @@ function SciMLBase.__solve(prob::NonlinearProblem, alg::SimpleTrustRegion, args.
121147 fx, ∇f = value_and_jacobian (autodiff, prob. f, fx, x, prob. p, jac_cache; J)
122148
123149 # Update the trust region radius.
124- (r > η₃) && (norm (δ) ≈ Δ) && (Δ = min (t₂ * Δ, Δₘₐₓ))
150+ if ! _unwrap_val (alg. nlsolve_update_rule) && r > η₃
151+ Δ = min (t₂ * Δ, Δₘₐₓ)
152+ end
125153 fₖ = fₖ₊₁
126154
127155 @bb H = transpose (∇f) × ∇f
128156 @bb g = transpose (∇f) × vec (fx)
129157 end
158+
159+ if _unwrap_val (alg. nlsolve_update_rule)
160+ if r > η₃
161+ Δ = t₂ * norm (δ)
162+ elseif r > 0.5
163+ Δ = max (Δ, t₂ * norm (δ))
164+ end
165+ end
130166 end
131167
132168 return build_solution (prob, alg, x, fx; retcode = ReturnCode. MaxIters)
0 commit comments