Adress points from code review.

Author: kellertuer

Changelog.md

@@ -20,14 +20,15 @@ We also unified a few of the internal solver state constructors.
 * New stopping criteria: `StopWhenRelativeAPosterioriCostChangeLessOrEqual` and `StopWhenProjectedNegativeGradientNormLess`. (#554).
 * `HagerZhangLinesearch` stepsize, a state-of-the-art line search for smooth objectives with cubic interpolation and adaptive Wolfe condition checking. (#554)
 * Stopping criteria can now be initialized using `initialize_stepsize!`, similar to solvers. (#554)
+* The `ConjugateResidualState` now has a `warm_start=` option when used multiple times, for example in every iteration as a subsolver, to reuse the last state from the previous run.
 
 ### Changed
 
 * In the [Riemannian Levenberg Marquardt algorithm](https://manoptjl.org/stable/solvers/LevenbergMarquardt/)t the `η` parameter has been renamed to `candidate_acceptance_threshold`, `β` to `damping_increase_factor` and `β_reduction` to `damping_reduction_factor`. (#617)
 * the constructor for the [Levenberg-Marquardt state](https://manoptjl.org/stable/solvers/LevenbergMarquardt/#Manopt.LevenbergMarquardtState) has been unified with the remaining states, to take the `sub_problem` and `sub_state` arguments as second and third positional arguments, respectively. (#617)
 * the keyword `initial_jacobian_f` within `LevenbergMarquardt` is unified in naming to the residual values vector and called `initial_jacobian_matrices`. If you call `LevenbergMarquardt` with a single vector component, also a single matrix is allowed. (#617)
 * an internal field of the solver state of Levenberg-Marqwuardt was called `jacobian_f` the same as the functions whose result it meant to cache if applicable. To distinguish both, the field is now called `jacobian_matrices`. (#617)
-* title of "How to define the cost in the embedding" tutorial (#615)
+* the `max_stepsize(M)` on the [`SymmetricPositiveDefinite`](https://juliamanifolds.github.io/Manifolds.jl/stable/manifolds/symmetricpositivedefinite/) manifold was changed from returning `Inf`, which is the mathematical maximal stepsize to returning the square root of the maximum (floating point) value to avoid numerical instabilities.
 
 ### Fixed
 

docs/src/notation.md

@@ -9,6 +9,6 @@ with the following additional parts.
 | ``\mathbf{A}`` | a matrix | ``\mathbf{B}, \mathbf{C}`` | |
 | ``\mathcal{A}, \mathcal{L}`` | linear operators, usually on tangent spaces | | |
 | ``\mathcal{C}`` | a (geodesically) convex set | | |
-| ``\operatorname*{arg\,min}_v`` | argument/variable ``v`` of a function ``f`` where a local or global minimum is attained | |
+| ``\operatorname*{arg\,min}_v f(v)`` | argument/variable ``v`` of a function ``f`` where a local or global minimum is attained | |
 | ``k`` | the current iterate | ``i`` | the goal is to unify this to `k` |
 | ``∇`` | The [Levi-Cevita connection](https://en.wikipedia.org/wiki/Levi-Civita_connection) | | |

docs/src/plans/objectives/constrained.md

@@ -9,8 +9,6 @@ ConstrainedManifoldObjective
 ManifoldConstrainedSetObjective
 ```
 
-It might be beneficial to use the adapted problem to specify different ranges for the gradients of the constraints
-
 ```@docs
 ConstrainedManoptProblem
 ```

src/plans/robustifiers.jl

@@ -30,7 +30,7 @@ Evaluate the robustifier ``ρ`` and its first two derivatives at `x`.
 A tuple `(a, b, c)` with
 * `a = ρ(x)`
 * `b = ρ'(x)`
-* `c = ρ''(x)` (might be `missing` if not defined)
+* `c = ρ''(x)`
 """
 get_robustifier_values(ρ::AbstractRobustifierFunction, x::Real)
 
@@ -52,17 +52,14 @@ and second derivative ``ρ'`` and ``ρ''``, respectively.
 
 Generate a `RobustifierFunction` given the function `ρ` and its first and second derivative.
 """
-struct RobustifierFunction{F <: Function, G <: Function, H <: Union{Function, Missing}} <: AbstractRobustifierFunction
+struct RobustifierFunction{F <: Function, G <: Function, H <: Function} <: AbstractRobustifierFunction
     ρ::F
     ρ_prime::G
     ρ_double_prime::H
 end
 
 function get_robustifier_values(ρf::RobustifierFunction, x::Real)
-    a = ρf.ρ(x)
-    b = ρf.ρ_prime(x)
-    c = ismissing(ρf.ρ_double_prime) ? missing : ρf.ρ_double_prime(x)
-    return (a, b, c)
+    return (ρf.ρ(x), ρf.ρ_prime(x), ρf.ρ_double_prime(x))
 end
 
 #
@@ -204,7 +201,7 @@ function get_robustifier_values(
     (a1, b1, c1) = get_robustifier_values(crf.ρ1, a2)
     a = a1
     b = b1 * b2
-    c = (ismissing(c1) || ismissing(c2)) ? missing : c1 * b2^2 + b1 * c2
+    c = c1 * b2^2 + b1 * c2
     return (a, b, c)
 end
 
@@ -314,7 +311,7 @@ function get_robustifier_values(srf::ScaledRobustifierFunction, x::Real)
     (a, b, c) = get_robustifier_values(srf.robustifier, x / s2)
     a_scaled = s2 * a
     b_scaled = b
-    c_scaled = ismissing(c) ? missing : c / s2
+    c_scaled = c / s2
     return (a_scaled, b_scaled, c_scaled)
 end
 

src/solvers/LevenbergMarquardt.jl

@@ -77,7 +77,7 @@ as well as in general using the model improvement parameter ``m_k`` in several p
 * `initial_jacobian_matrices`: a cache for the evaluated Jacobians (currently only used if `use_unified_basis = true`, then initialised to a vector of jacobian matrices, otherwise ignored)
 $(_kwargs(:retraction_method))
 * `scaling_threshold = 1.0e-6`:         a threshold `ε` to bound the scaling parameter `α` in the robust case away from `1`, see [`get_LevenbergMarquardt_scaling`](@ref)
-* `scaling_mode = :Default`:            specify the scaling stabilization mode, see [`get_LevenbergMarquardt_scaling`](@ref)
+* `scaling_mode = :Strict`:            specify the scaling stabilization mode, see [`get_LevenbergMarquardt_scaling`](@ref)
 $(_kwargs(:stopping_criterion; default = "`[`StopAfterIteration`](@ref)`(500)`$(_sc(:Any))[`StopWhenGradientNormLess`](@ref)`(1.0e-12)`$(_sc(:Any))[`StopWhenStepsizeLess`](@ref)`(1.0e-12)"))
 * `sub_objective`                      : specify the objective for the surrogate sub problem to solver in every iteration.
   This is set depending on the `use_unified_basis`
@@ -238,7 +238,7 @@ function LevenbergMarquardt!(
             fill(nothing, length(get_objective(nlso).objective))
         end,
         scaling_threshold::Real = 1.0e-6,
-        scaling_mode::Symbol = :Default,
+        scaling_mode::Symbol = :Strict,
         minimum_acceptable_model_improvement::Real = eps(number_eltype(p)),
         sub_objective = construct_lm_subobjective(use_unified_basis, nlso, damping_term_min, scaling_threshold, scaling_mode, initial_residual_values, initial_jacobian_matrices),
         sub_problem = DefaultManoptProblem(TangentSpace(M, p), sub_objective),