lin ineq constraint value cache fix

src/SimpleMOP/simple_caches.jl

@@ -18,11 +18,13 @@ Base.@kwdef struct SimpleValueCache{
 
 	Ax_min_b :: Vector{F}
 	Ex_min_c :: Vector{F}
-
+
+    #=
 	"Reference to maximum constraint violation."
     theta_ref :: Base.RefValue{F}
 	"Reference to maximum function value."
     phi_ref :: Base.RefValue{F}
+    =#
 end
 
 struct SimpleSurrogateValueCache{
@@ -83,14 +85,15 @@ function init_simple_value_caches_for_mop(T, nx, nfx, nhx, ngx, nE, nA)
     Ax_min_b = similar(Ax)
     Ex_min_c = similar(Ex)
 
+    #=
     ## constraint violation and filter value
-    theta_ref = Ref(zero(T))
-    phi_ref = Ref(zero(T))
-
+    theta_ref = Ref(T(NaN))
+    phi_ref = Ref(T(NaN))
+    =#
     return SimpleValueCache(; 
         ξ, x, fx, hx, gx, 
         Ax, Ex, Ax_min_b, Ex_min_c, 
-        theta_ref, phi_ref
+    #    theta_ref, phi_ref
     )
 end
 

src/mop.jl

@@ -180,7 +180,7 @@ function eval_mop!(mop_cache, mop)
         cached_Ax(mop_cache), 
         mop, ξ
     )
-    return nothing
+    return θ, Φ
 end
 
 "Evaluate `mop` at unscaled site `ξ` and modify result arrays in place."

src/multi_filters.jl

@@ -377,7 +377,9 @@ is_converged(sstructs::SolutionStructs)=is_converged(sstructs.status_ref)
 @forward SolutionStructs.vals cached_Ax(sols::SolutionStructs)
 @forward SolutionStructs.vals cached_Ex(sols::SolutionStructs)
 @forward SolutionStructs.vals cached_Ex_min_c(sols::SolutionStructs)
+@forward SolutionStructs.vals cached_Ax_min_b(sols::SolutionStructs)
 @forward SolutionStructs.vals cached_theta(sols::SolutionStructs)
+@forward SolutionStructs.vals cached_Phi(sols::SolutionStructs)
 dim_theta(::SolutionStructs) = 1
 
 function copy_solution_structs(sstructs)

src/set_algo.jl

@@ -647,6 +647,8 @@ function optimize_many(
         iteration_status, iteration_scalars,
     ) = optimizer_caches
 
+    # for the moment, only use minimum radius stopping (and max_iter)
+    # TODO re-enable other criteria
     @unpack stop_delta_min = algo_opts
     stop_crits = MinimumRadiusStopping(; delta_min = stop_delta_min)
     @reset optimizer_caches.stop_crits = stop_crits
@@ -888,7 +890,12 @@ function _test_trial_point!(
             $(indent_str(indent))- mx - mxs = $(pretty_row_vec(diff_fx_mod)),
             $(indent_str(indent))- theta(xs)= $(cached_theta(vals_tmp)),
             $(indent_str(indent))- κθ^ψ     = $(kappa_theta * θx^psi_theta)"""
-
+
+
+        ## `is_pos` is an index vector for the positive entries of the difference
+        ## vector `del_mx`.
+        ## We use it in case of partial descent.
+        ## `Δmx` is the minimum surrogate objective reduction.
         is_pos = zeros(Bool, length(diff_fx_mod))
         Δmx = Inf
         for (l, del_mx) in enumerate(diff_fx_mod)
@@ -906,31 +913,40 @@ function _test_trial_point!(
             is_good_trial_point = false
             delta_child = delta_new = gamma_shrink_much * delta
         else
-            #=
-            ## if m(x) ⪨ m(x+s) + κ⋅θ(x)^ψ
-            mxs_θ =  mxs .+(kappa_theta * θx^psi_theta)
-            if vector_dominates(mx[is_pos], mxs_θ[is_pos])
-            =#
-            if !isinf(Δmx) && Δmx < kappa_theta * θx^psi_theta
+            ## The single-objective f-step test is 
+            ## m(x) - m(x+s) ≥ κ θ(x)^ψ
+            ## ⇔
+            ## m(x + s) + κ θ(x)^ψ ≤ m(x)
+            ## In the multi-objective case, let's try
+            ## m(x+s) + κθ(x)^ψ ≼ m(x).
+            ## If this test fails, we deem the constraint violation too high and
+            ## have a θ-step.
+            ## For sake of simplicity, we test 
+            ## any( m(x+s) + κθ(x)^ψ .> m(x) )
+            ## If that is true, then we have θ-step.
+            ##
+            ## In case of partial descent, we have to respect the index set `is_pos`.
+            theta_step_test_offset = kappa_theta * θx^psi_theta
+            if isinf(Δmx) || any( mxs[is_pos] .+ theta_step_test_offset .> mx[is_pos] )
                 ## (theta step)
                 ## trial point is accepted    
                 is_good_trial_point = true
+                ## do not change radius # TODO could be smarter
                 delta_child = delta_new = delta
                 ## sol is added to filter (sync population!)
                 ## (adding sol to filter makes it incompatible with filter, 
-                ## it is removed from population
+                ##   it should be/is (?) removed from population)
                 augment_filter = true
             else
-                offset = Δmx * nu_accept
-                if !isinf(Δmx) && trial_point_population_acceptable(fxs, θxs, population, offset)
+                offset = zero(fxs)  # TODO cache
+                offset[is_pos] .= diff_fx_mod[is_pos] .* nu_accept
+                if trial_point_population_acceptable(fxs, θxs, population, offset)
                     is_good_trial_point = true
 
                     delta_child = delta_new = gamma_shrink * delta
-                    rho_success = minimum( diff_fx[is_pos] ./ diff_fx_mod[is_pos] )
-                    if all(is_pos) && rho_success >= nu_success
-                        if delta < delta_max 
-                            delta_child = min(delta_max, gamma_grow * delta)
-                        end
+                    rho_success = minimum( diff_fx[is_pos] ./ diff_fx_mod[is_pos]; init=-Inf )
+                    if rho_success >= nu_success
+                        delta_child = min(delta_max, gamma_grow * delta)
                     end
                 else
                     ## (innacceptable)
@@ -956,11 +972,11 @@ function cosine_similarity(p1, p2)
 end
 
 function trial_point_population_acceptable(fxs, θxs, population, offset=0)
-    fxs = fxs .+ offset     # TODO cache
+    fxs_shifted = fxs .+ offset     # TODO cache
     for sol in nondominated_elements(population)
         if augmented_dominates(
             cached_theta(sol), cached_fx(sol),
-            θxs, fxs, 
+            θxs, fxs_shifted, 
         )
             return false
         end

src/value_caches.jl

@@ -37,17 +37,20 @@ dim_lin_ineq_constraints(c::AbstractValueCache)=length(cached_Ax(c))
 
 struct WrappedMOPCache{F, cacheType<:AbstractMOPCache{F}}<:AbstractMOPCache{F}
 	cache :: cacheType
-	xhash :: Base.RefValue{UInt64}
+	hash_theta :: Base.RefValue{UInt64}
+	hash_Phi :: Base.RefValue{UInt64}
 	theta_ref :: Base.RefValue{F}
 	Phi_ref :: Base.RefValue{F}
 end
 
 function WrappedMOPCache(cache::AbstractMOPCache{F}) where F
+    NaNF = F(NaN)
     return WrappedMOPCache(
         cache,
-        Ref(zero(UInt64)),
-        Ref(F(NaN)),
-        Ref(F(NaN)),
+        Ref(hash(NaNF)),
+        Ref(hash(NaNF)),
+        Ref(NaNF),
+        Ref(NaNF),
     )
 end
 
@@ -68,35 +71,33 @@ end
 @forward WrappedMOPCache.cache dim_lin_eq_constraints(wcache::WrappedMOPCache)
 @forward WrappedMOPCache.cache dim_lin_ineq_constraints(wcache::WrappedMOPCache)
 
-function cached_theta(wcache::WrappedMOPCache, step_down=true)
+function cached_theta(wcache::WrappedMOPCache)
     cache = wcache.cache
     x = cached_x(cache)
-    xhash = wcache.xhash[]
-    inhash = hash(x)
-    if xhash != inhash
+    hash_old = wcache.hash_theta[]
+    hash_x = hash(x)
+    if hash_old != hash_x
         theta = constraint_violation(
             cached_hx(cache),
             cached_gx(cache),
             cached_Ex_min_c(cache),
             cached_Ax_min_b(cache)
         )
-        step_down && cached_Phi(wcache, false)
         wcache.theta_ref[] = theta
-        wcache.xhash[] = inhash
+        wcache.hash_theta[] = hash_x
     end
     return wcache.theta_ref[]
 end
 
-function cached_Phi(wcache::WrappedMOPCache, step_down=true)
+function cached_Phi(wcache::WrappedMOPCache)
     cache = wcache.cache
     x = cached_x(cache)
-    xhash = wcache.xhash[]
-    inhash = hash(x)
-    if xhash != inhash
+    hash_old = wcache.hash_Phi[]
+    hash_x = hash(x)
+    if hash_old != hash_x
         Phi = maximum(cached_fx(cache))
-        step_down && cached_theta(wcache, false)
         wcache.Phi_ref[] = Phi
-        wcache.xhash[] = inhash
+        wcache.hash_Phi[] = hash_x
     end
     return wcache.Phi_ref[]
 end

test/ndset.jl

@@ -4,7 +4,7 @@ import Compromise: @unpack
 includet("TestProblems.jl")
 TP = TestProblems
 
-tp = TP._test_problem(Val(3), 2)
+tp = TP._test_problem(Val(6), 2)
 _mop = TP.to_mutable_mop(tp; max_func_calls=10000)
 X = [3.42e-01; 9.23e-01;;]
 X = _mop.lb .+ (_mop.ub .- _mop.lb) .* rand(2, 50)
@@ -14,7 +14,7 @@ algo_opts = C.AlgorithmOptions(;
     max_iter=100, 
     stop_delta_min=1e-5,
     step_config = C.SteepestDescentConfig(;
-        backtracking_mode = Val(:all)
+        backtracking_mode = Val(:any)
     )
 )