Skip to main content

oxilean_std/convex_analysis/
functions.rs

1//! Auto-generated module
2//!
3//! 🤖 Generated with [SplitRS](https://github.com/cool-japan/splitrs)
4
5use oxilean_kernel::{BinderInfo, Declaration, Environment, Expr, Level, Name};
6
7use super::types::{
8    ADMMData, AdmmSolver, AlternatingProjectionSolver, ConvexCone, ConvexConjugate, ConvexFunction,
9    ConvexProblemClass, ConvexProgram, Epigraph, ExtragradientMethod, FenchelConjugate,
10    FenchelConjugateEvaluator, FenchelDualityPair, FunctionClass, Hyperplane, LagrangianDuality,
11    MirrorDescent, ProxConfig, ProximalOperator, ProximalOperatorNew, ProximalPointAlgorithm,
12    RecessionCone, SeparatingHyperplaneFinder, StepSchedule, Subdifferential, SubgradientMethod,
13    SublevelSet, VariationalInequality,
14};
15
16pub fn app(f: Expr, a: Expr) -> Expr {
17    Expr::App(Box::new(f), Box::new(a))
18}
19pub fn app2(f: Expr, a: Expr, b: Expr) -> Expr {
20    app(app(f, a), b)
21}
22pub fn app3(f: Expr, a: Expr, b: Expr, c: Expr) -> Expr {
23    app(app2(f, a, b), c)
24}
25pub fn cst(s: &str) -> Expr {
26    Expr::Const(Name::str(s), vec![])
27}
28pub fn prop() -> Expr {
29    Expr::Sort(Level::zero())
30}
31pub fn type0() -> Expr {
32    Expr::Sort(Level::succ(Level::zero()))
33}
34pub fn pi(bi: BinderInfo, name: &str, dom: Expr, body: Expr) -> Expr {
35    Expr::Pi(bi, Name::str(name), Box::new(dom), Box::new(body))
36}
37pub fn arrow(a: Expr, b: Expr) -> Expr {
38    pi(BinderInfo::Default, "_", a, b)
39}
40pub fn bvar(n: u32) -> Expr {
41    Expr::BVar(n)
42}
43pub fn real_ty() -> Expr {
44    cst("Real")
45}
46pub fn nat_ty() -> Expr {
47    cst("Nat")
48}
49pub fn list_ty(elem: Expr) -> Expr {
50    app(cst("List"), elem)
51}
52pub fn bool_ty() -> Expr {
53    cst("Bool")
54}
55pub fn fn_ty(dom: Expr, cod: Expr) -> Expr {
56    arrow(dom, cod)
57}
58pub fn vec_ty() -> Expr {
59    list_ty(real_ty())
60}
61pub fn mat_ty() -> Expr {
62    list_ty(list_ty(real_ty()))
63}
64pub fn extended_real_ty() -> Expr {
65    cst("ExtendedReal")
66}
67/// `IsConvexSet : (List Real -> Prop) -> Prop`
68/// A set C ⊆ ℝ^n given as characteristic predicate is convex.
69pub fn is_convex_set_ty() -> Expr {
70    arrow(fn_ty(vec_ty(), prop()), prop())
71}
72/// `IsConvexFunction : (List Real -> Real) -> Prop`
73/// f is convex: f(λx + (1-λ)y) ≤ λf(x) + (1-λ)f(y) for all λ ∈ [0,1].
74pub fn is_convex_function_ty() -> Expr {
75    arrow(fn_ty(vec_ty(), real_ty()), prop())
76}
77/// `IsStrictlyConvexFunction : (List Real -> Real) -> Prop`
78/// f is strictly convex: strict inequality for λ ∈ (0,1) and x ≠ y.
79pub fn is_strictly_convex_function_ty() -> Expr {
80    arrow(fn_ty(vec_ty(), real_ty()), prop())
81}
82/// `IsStronglyConvexFunction : (List Real -> Real) -> Real -> Prop`
83/// f is strongly convex with modulus m: f(λx+(1-λ)y) ≤ λf(x)+(1-λ)f(y) - (m/2)λ(1-λ)‖x-y‖².
84pub fn is_strongly_convex_function_ty() -> Expr {
85    arrow(fn_ty(vec_ty(), real_ty()), arrow(real_ty(), prop()))
86}
87/// `Epigraph : (List Real -> Real) -> (List Real -> Prop)`
88/// epi(f) = {(x, t) | f(x) ≤ t}, encoded as a predicate on ℝ^{n+1}.
89pub fn epigraph_ty() -> Expr {
90    arrow(fn_ty(vec_ty(), real_ty()), fn_ty(vec_ty(), prop()))
91}
92/// `LevelSet : (List Real -> Real) -> Real -> (List Real -> Prop)`
93/// lev_α(f) = {x | f(x) ≤ α}.
94pub fn level_set_ty() -> Expr {
95    arrow(
96        fn_ty(vec_ty(), real_ty()),
97        arrow(real_ty(), fn_ty(vec_ty(), prop())),
98    )
99}
100/// `ClosureFunction : (List Real -> Real) -> (List Real -> Real)`
101/// cl(f): lower semi-continuous closure of f.
102pub fn closure_function_ty() -> Expr {
103    arrow(fn_ty(vec_ty(), real_ty()), fn_ty(vec_ty(), real_ty()))
104}
105/// `IsLowerSemicontinuous : (List Real -> Real) -> Prop`
106/// f is lsc iff its epigraph is closed.
107pub fn is_lsc_ty() -> Expr {
108    arrow(fn_ty(vec_ty(), real_ty()), prop())
109}
110/// `IsProperConvex : (List Real -> Real) -> Prop`
111/// f is proper: never −∞ and not identically +∞.
112pub fn is_proper_convex_ty() -> Expr {
113    arrow(fn_ty(vec_ty(), real_ty()), prop())
114}
115/// `SupportingHyperplane : (List Real -> Prop) -> List Real -> Prop`
116/// At boundary point x₀ of convex set C, ∃ nonzero c with c·x ≤ c·x₀ for all x ∈ C.
117pub fn supporting_hyperplane_ty() -> Expr {
118    arrow(fn_ty(vec_ty(), prop()), arrow(vec_ty(), prop()))
119}
120/// `SeparatingHyperplane : (List Real -> Prop) -> (List Real -> Prop) -> Prop`
121/// Two disjoint convex sets can be separated by a hyperplane.
122pub fn separating_hyperplane_ty() -> Expr {
123    let set_ty = fn_ty(vec_ty(), prop());
124    arrow(set_ty.clone(), arrow(set_ty, prop()))
125}
126/// `StrictSeparation : (List Real -> Prop) -> (List Real -> Prop) -> Prop`
127/// Compact convex C and closed convex D, disjoint, can be strictly separated.
128pub fn strict_separation_ty() -> Expr {
129    let set_ty = fn_ty(vec_ty(), prop());
130    arrow(set_ty.clone(), arrow(set_ty, prop()))
131}
132/// `SupportFunction : (List Real -> Prop) -> List Real -> Real`
133/// σ_C(y) = sup_{x ∈ C} ⟨y, x⟩.
134pub fn support_function_ty() -> Expr {
135    arrow(fn_ty(vec_ty(), prop()), fn_ty(vec_ty(), real_ty()))
136}
137/// `GaugeFunction : (List Real -> Prop) -> List Real -> Real`
138/// γ_C(x) = inf{t ≥ 0 | x ∈ t·C}.
139pub fn gauge_function_ty() -> Expr {
140    arrow(fn_ty(vec_ty(), prop()), fn_ty(vec_ty(), real_ty()))
141}
142/// `FenchelConjugate : (List Real -> Real) -> (List Real -> Real)`
143/// f*(y) = sup_x {⟨y,x⟩ − f(x)}.
144pub fn fenchel_conjugate_ty() -> Expr {
145    let rn_r = fn_ty(vec_ty(), real_ty());
146    arrow(rn_r.clone(), rn_r)
147}
148/// `Biconjugate : (List Real -> Real) -> (List Real -> Real)`
149/// f** = cl(conv(f)): the lsc convex hull of f.
150pub fn biconjugate_ty() -> Expr {
151    let rn_r = fn_ty(vec_ty(), real_ty());
152    arrow(rn_r.clone(), rn_r)
153}
154/// `FenchelYoungInequality : (List Real -> Real) -> Prop`
155/// For any f and its conjugate f*: ⟨x, y⟩ ≤ f(x) + f*(y).
156pub fn fenchel_young_inequality_ty() -> Expr {
157    arrow(fn_ty(vec_ty(), real_ty()), prop())
158}
159/// `FenchelYoungEquality : (List Real -> Real) -> Prop`
160/// ⟨x, y⟩ = f(x) + f*(y) iff y ∈ ∂f(x) (subgradient condition).
161pub fn fenchel_young_equality_ty() -> Expr {
162    arrow(fn_ty(vec_ty(), real_ty()), prop())
163}
164/// `LegendreFenchelTransform : (List Real -> Real) -> (List Real -> Real)`
165/// Alias for FenchelConjugate for convex lsc functions.
166pub fn legendre_fenchel_transform_ty() -> Expr {
167    fenchel_conjugate_ty()
168}
169/// `ConjugateOfSum : Prop`
170/// (f + g)*(y) = (f* □ g*)(y) where □ is inf-convolution (under qualification).
171pub fn conjugate_of_sum_ty() -> Expr {
172    prop()
173}
174/// `MorozovIdentity : (List Real -> Real) -> Prop`
175/// Moreau decomposition: prox_f(x) + prox_{f*}(x) = x.
176pub fn moreau_identity_ty() -> Expr {
177    arrow(fn_ty(vec_ty(), real_ty()), prop())
178}
179/// `Subgradient : (List Real -> Real) -> List Real -> List Real -> Prop`
180/// g ∈ ∂f(x) iff f(y) ≥ f(x) + ⟨g, y - x⟩ for all y.
181pub fn subgradient_ty() -> Expr {
182    arrow(
183        fn_ty(vec_ty(), real_ty()),
184        arrow(vec_ty(), arrow(vec_ty(), prop())),
185    )
186}
187/// `Subdifferential : (List Real -> Real) -> List Real -> (List Real -> Prop)`
188/// ∂f(x) = {g | g is a subgradient of f at x}.
189pub fn subdifferential_ty() -> Expr {
190    arrow(
191        fn_ty(vec_ty(), real_ty()),
192        arrow(vec_ty(), fn_ty(vec_ty(), prop())),
193    )
194}
195/// `NormalCone : (List Real -> Prop) -> List Real -> (List Real -> Prop)`
196/// N_C(x) = {v | ⟨v, y-x⟩ ≤ 0 for all y ∈ C}.
197pub fn normal_cone_ty() -> Expr {
198    arrow(
199        fn_ty(vec_ty(), prop()),
200        arrow(vec_ty(), fn_ty(vec_ty(), prop())),
201    )
202}
203/// `TangentCone : (List Real -> Prop) -> List Real -> (List Real -> Prop)`
204/// T_C(x) = cl({d | ∃ t_k ↘ 0, x + t_k d ∈ C}).
205pub fn tangent_cone_ty() -> Expr {
206    arrow(
207        fn_ty(vec_ty(), prop()),
208        arrow(vec_ty(), fn_ty(vec_ty(), prop())),
209    )
210}
211/// `OptimalityConditionConvex : (List Real -> Real) -> (List Real -> Prop) -> Prop`
212/// x* minimizes f over C iff 0 ∈ ∂f(x*) + N_C(x*).
213pub fn optimality_condition_convex_ty() -> Expr {
214    arrow(
215        fn_ty(vec_ty(), real_ty()),
216        arrow(fn_ty(vec_ty(), prop()), prop()),
217    )
218}
219/// `SubdiffOfSum : Prop`
220/// Under constraint qualification: ∂(f+g)(x) ⊇ ∂f(x) + ∂g(x) (Moreau-Rockafellar).
221pub fn subdiff_of_sum_ty() -> Expr {
222    prop()
223}
224/// `SubdiffOfComposition : Prop`
225/// Chain rule for subdifferentials of composed convex functions.
226pub fn subdiff_of_composition_ty() -> Expr {
227    prop()
228}
229/// `RecessionCone : (List Real -> Prop) -> (List Real -> Prop)`
230/// rec(C) = {d | x + td ∈ C for all t ≥ 0, for some x ∈ C}.
231pub fn recession_cone_ty() -> Expr {
232    let set_ty = fn_ty(vec_ty(), prop());
233    arrow(set_ty.clone(), set_ty)
234}
235/// `RecessionFunction : (List Real -> Real) -> (List Real -> Real)`
236/// rec(f)(d) = lim_{t→∞} f(x + td)/t (horizon function of f).
237pub fn recession_function_ty() -> Expr {
238    let rn_r = fn_ty(vec_ty(), real_ty());
239    arrow(rn_r.clone(), rn_r)
240}
241/// `IsCoercive : (List Real -> Real) -> Prop`
242/// f is coercive: f(x) → +∞ as ‖x‖ → ∞.
243pub fn is_coercive_ty() -> Expr {
244    arrow(fn_ty(vec_ty(), real_ty()), prop())
245}
246/// `IsLevelBounded : (List Real -> Real) -> Prop`
247/// All level sets {x | f(x) ≤ α} are bounded.
248pub fn is_level_bounded_ty() -> Expr {
249    arrow(fn_ty(vec_ty(), real_ty()), prop())
250}
251/// `ExistenceOfMinimum : (List Real -> Real) -> (List Real -> Prop) -> Prop`
252/// Weierstrass theorem: lsc proper coercive f on closed set has a minimizer.
253pub fn existence_of_minimum_ty() -> Expr {
254    arrow(
255        fn_ty(vec_ty(), real_ty()),
256        arrow(fn_ty(vec_ty(), prop()), prop()),
257    )
258}
259/// `ProximalOperator : (List Real -> Real) -> Real -> List Real -> List Real`
260/// prox_{λf}(v) = argmin_x {f(x) + (1/2λ)‖x-v‖²}.
261pub fn proximal_operator_ty() -> Expr {
262    arrow(
263        fn_ty(vec_ty(), real_ty()),
264        arrow(real_ty(), fn_ty(vec_ty(), vec_ty())),
265    )
266}
267/// `MoreauEnvelope : (List Real -> Real) -> Real -> List Real -> Real`
268/// e_{λf}(x) = inf_y {f(y) + (1/2λ)‖x-y‖²}.
269pub fn moreau_envelope_ty() -> Expr {
270    arrow(
271        fn_ty(vec_ty(), real_ty()),
272        arrow(real_ty(), fn_ty(vec_ty(), real_ty())),
273    )
274}
275/// `InfConvolution : (List Real -> Real) -> (List Real -> Real) -> List Real -> Real`
276/// (f □ g)(x) = inf_{y} {f(y) + g(x - y)}.
277pub fn inf_convolution_ty() -> Expr {
278    let rn_r = fn_ty(vec_ty(), real_ty());
279    arrow(rn_r.clone(), arrow(rn_r, fn_ty(vec_ty(), real_ty())))
280}
281/// `ProxFirmlyNonexpansive : (List Real -> Real) -> Prop`
282/// prox_{λf} is firmly nonexpansive: ‖prox x - prox y‖² ≤ ⟨prox x - prox y, x - y⟩.
283pub fn prox_firmly_nonexpansive_ty() -> Expr {
284    arrow(fn_ty(vec_ty(), real_ty()), prop())
285}
286/// `MoreauEnvelopeSmooth : (List Real -> Real) -> Prop`
287/// e_{λf} is 1/λ-smooth and ∇e_{λf}(x) = (1/λ)(x - prox_{λf}(x)).
288pub fn moreau_envelope_smooth_ty() -> Expr {
289    arrow(fn_ty(vec_ty(), real_ty()), prop())
290}
291/// `MoreauEnvelopeConvergence : (List Real -> Real) -> Prop`
292/// e_{λf}(x) ↗ f(x) as λ ↘ 0 (pointwise approximation).
293pub fn moreau_envelope_convergence_ty() -> Expr {
294    arrow(fn_ty(vec_ty(), real_ty()), prop())
295}
296/// `BregmanDivergence : (List Real -> Real) -> List Real -> List Real -> Real`
297/// D_f(x, y) = f(x) - f(y) - ⟨∇f(y), x - y⟩ for differentiable f.
298pub fn bregman_divergence_ty() -> Expr {
299    arrow(
300        fn_ty(vec_ty(), real_ty()),
301        fn_ty(vec_ty(), fn_ty(vec_ty(), real_ty())),
302    )
303}
304/// `BregmanProjection : (List Real -> Real) -> (List Real -> Prop) -> List Real -> List Real`
305/// argmin_{y ∈ C} D_f(y, x) — Bregman projection of x onto C.
306pub fn bregman_projection_ty() -> Expr {
307    arrow(
308        fn_ty(vec_ty(), real_ty()),
309        arrow(fn_ty(vec_ty(), prop()), fn_ty(vec_ty(), vec_ty())),
310    )
311}
312/// `BregmanThreePointIdentity : (List Real -> Real) -> Prop`
313/// D_f(x, z) = D_f(x, y) + D_f(y, z) + ⟨∇f(y) - ∇f(z), x - y⟩.
314pub fn bregman_three_point_identity_ty() -> Expr {
315    arrow(fn_ty(vec_ty(), real_ty()), prop())
316}
317/// `MirrorDescentStep : (List Real -> Real) -> (List Real -> Real) -> Real -> List Real -> List Real`
318/// Mirror descent update: x_{k+1} = argmin_{y} {⟨g_k, y⟩ + (1/η) D_f(y, x_k)}.
319pub fn mirror_descent_step_ty() -> Expr {
320    let rn_r = fn_ty(vec_ty(), real_ty());
321    arrow(
322        rn_r.clone(),
323        arrow(rn_r, arrow(real_ty(), fn_ty(vec_ty(), vec_ty()))),
324    )
325}
326/// `ExtremePoint : (List Real -> Prop) -> List Real -> Prop`
327/// x is an extreme point of C: x ∉ (y, z) strictly for any y, z ∈ C with y ≠ z.
328pub fn extreme_point_ty() -> Expr {
329    arrow(fn_ty(vec_ty(), prop()), arrow(vec_ty(), prop()))
330}
331/// `CaratheodoryTheorem : (List Real -> Prop) -> Nat -> Prop`
332/// Every point in conv(C) ⊆ ℝ^n can be written as a convex combination of ≤ n+1 points of C.
333pub fn caratheodory_theorem_ty() -> Expr {
334    arrow(fn_ty(vec_ty(), prop()), arrow(nat_ty(), prop()))
335}
336/// `KreinMilmanTheorem : (List Real -> Prop) -> Prop`
337/// Every compact convex set is the closed convex hull of its extreme points.
338pub fn krein_milman_theorem_ty() -> Expr {
339    arrow(fn_ty(vec_ty(), prop()), prop())
340}
341/// `ConvexHull : (List Real -> Prop) -> (List Real -> Prop)`
342/// conv(C) = smallest convex set containing C.
343pub fn convex_hull_ty() -> Expr {
344    let set_ty = fn_ty(vec_ty(), prop());
345    arrow(set_ty.clone(), set_ty)
346}
347/// `ClosedConvexHull : (List Real -> Prop) -> (List Real -> Prop)`
348/// cl(conv(C)).
349pub fn closed_convex_hull_ty() -> Expr {
350    let set_ty = fn_ty(vec_ty(), prop());
351    arrow(set_ty.clone(), set_ty)
352}
353/// `FarkasLemma : Prop`
354/// Ax = b, x ≥ 0 has a solution iff for all y with A^Ty ≥ 0 we have b^Ty ≥ 0.
355pub fn farkas_lemma_ty() -> Expr {
356    prop()
357}
358/// `FarkasLemmaMatrix : (List (List Real)) -> List Real -> Prop`
359/// Farkas lemma parameterised by matrix A and vector b.
360pub fn farkas_lemma_matrix_ty() -> Expr {
361    arrow(mat_ty(), arrow(vec_ty(), prop()))
362}
363/// `GordonAlternative : (List (List Real)) -> Prop`
364/// Gordon's theorem: either ∃x, Ax < 0 or ∃y ≥ 0, y ≠ 0, A^Ty = 0, but not both.
365pub fn gordon_alternative_ty() -> Expr {
366    arrow(mat_ty(), prop())
367}
368/// `TuckerAlternative : (List (List Real)) -> Prop`
369/// Tucker's theorem of the alternative for linear inequalities.
370pub fn tucker_alternative_ty() -> Expr {
371    arrow(mat_ty(), prop())
372}
373/// `StiemkePosanAlternative : Prop`
374/// Stiemke/Positivstellensatz alternative for polynomial systems.
375pub fn stiemke_alternative_ty() -> Expr {
376    prop()
377}
378/// `IsLipschitzGradient : (List Real -> Real) -> Real -> Prop`
379/// ‖∇f(x) - ∇f(y)‖ ≤ L‖x-y‖ for all x, y.
380pub fn is_lipschitz_gradient_ty() -> Expr {
381    arrow(fn_ty(vec_ty(), real_ty()), arrow(real_ty(), prop()))
382}
383/// `StrongConvexityQuadraticBound : (List Real -> Real) -> Real -> Prop`
384/// f(y) ≥ f(x) + ⟨∇f(x), y-x⟩ + (m/2)‖y-x‖² for all x, y.
385pub fn strong_convexity_quadratic_bound_ty() -> Expr {
386    arrow(fn_ty(vec_ty(), real_ty()), arrow(real_ty(), prop()))
387}
388/// `IsSmooth : (List Real -> Real) -> Real -> Prop`
389/// f is L-smooth: has L-Lipschitz gradient.
390pub fn is_smooth_ty() -> Expr {
391    arrow(fn_ty(vec_ty(), real_ty()), arrow(real_ty(), prop()))
392}
393/// `DescentLemma : (List Real -> Real) -> Real -> Prop`
394/// f(y) ≤ f(x) + ⟨∇f(x), y-x⟩ + (L/2)‖y-x‖² (upper quadratic bound for L-smooth f).
395pub fn descent_lemma_ty() -> Expr {
396    arrow(fn_ty(vec_ty(), real_ty()), arrow(real_ty(), prop()))
397}
398/// `ConvergenceRateGradient : (List Real -> Real) -> Real -> Prop`
399/// Gradient descent on L-smooth convex f converges at rate O(1/k).
400pub fn convergence_rate_gradient_ty() -> Expr {
401    arrow(fn_ty(vec_ty(), real_ty()), arrow(real_ty(), prop()))
402}
403/// `LinearConvergenceStrongConvex : (List Real -> Real) -> Real -> Real -> Prop`
404/// Gradient descent on m-strongly convex L-smooth f converges at linear rate (1 - m/L)^k.
405pub fn linear_convergence_strong_convex_ty() -> Expr {
406    arrow(
407        fn_ty(vec_ty(), real_ty()),
408        arrow(real_ty(), arrow(real_ty(), prop())),
409    )
410}
411/// `IsBarrierFunction : (List Real -> Prop) -> (List Real -> Real) -> Prop`
412/// B is a barrier for C: B(x) → +∞ as x approaches the boundary of C.
413pub fn is_barrier_function_ty() -> Expr {
414    arrow(
415        fn_ty(vec_ty(), prop()),
416        arrow(fn_ty(vec_ty(), real_ty()), prop()),
417    )
418}
419/// `IsSelfConcordant : (List Real -> Real) -> Prop`
420/// f is self-concordant: |∇³f(x)[d,d,d]| ≤ 2(∇²f(x)[d,d])^(3/2) for all d.
421pub fn is_self_concordant_ty() -> Expr {
422    arrow(fn_ty(vec_ty(), real_ty()), prop())
423}
424/// `LogBarrier : (List Real -> Prop) -> (List Real -> Real)`
425/// B(x) = -∑ log(constraints(x)) — standard log barrier for a polytope.
426pub fn log_barrier_ty() -> Expr {
427    arrow(fn_ty(vec_ty(), prop()), fn_ty(vec_ty(), real_ty()))
428}
429/// `CentralPath : (List Real -> Real) -> (List Real -> Prop) -> Real -> (List Real -> Prop)`
430/// x*(t) = argmin_x { t * f(x) + B(x) } — central path parameterised by t.
431pub fn central_path_ty() -> Expr {
432    arrow(
433        fn_ty(vec_ty(), real_ty()),
434        arrow(
435            fn_ty(vec_ty(), prop()),
436            arrow(real_ty(), fn_ty(vec_ty(), prop())),
437        ),
438    )
439}
440/// `InteriorPointConvergence : (List Real -> Real) -> (List Real -> Prop) -> Prop`
441/// Central path converges to a solution as t → ∞.
442pub fn interior_point_convergence_ty() -> Expr {
443    arrow(
444        fn_ty(vec_ty(), real_ty()),
445        arrow(fn_ty(vec_ty(), prop()), prop()),
446    )
447}
448/// `IsSOCConstraint : List Real -> Real -> Prop`
449/// ‖Ax + b‖ ≤ c^T x + d — second-order cone (Lorentz cone) constraint.
450pub fn is_soc_constraint_ty() -> Expr {
451    arrow(vec_ty(), arrow(real_ty(), prop()))
452}
453/// `LorentzCone : Nat -> (List Real -> Prop)`
454/// K_n = {(x, t) ∈ ℝ^{n+1} | ‖x‖ ≤ t}.
455pub fn lorentz_cone_ty() -> Expr {
456    arrow(nat_ty(), fn_ty(vec_ty(), prop()))
457}
458/// `PositiveSemidefiniteCone : Nat -> (List (List Real) -> Prop)`
459/// PSD cone: S_n^+ = {X ∈ S_n | X ⪰ 0}.
460pub fn psd_cone_ty() -> Expr {
461    arrow(nat_ty(), fn_ty(mat_ty(), prop()))
462}
463/// `SOCPDuality : Prop`
464/// Strong duality for second-order cone programming under Slater's condition.
465pub fn socp_duality_ty() -> Expr {
466    prop()
467}
468/// `SDPDuality : Prop`
469/// Strong duality for semidefinite programming under Slater's condition.
470pub fn sdp_duality_ty() -> Expr {
471    prop()
472}
473/// `LagrangianFunction : (List Real -> Real) -> (List (List Real -> Real)) -> List Real -> List Real -> Real`
474/// L(x, λ) = f(x) + ∑_i λ_i g_i(x) — Lagrangian of the primal problem.
475pub fn lagrangian_function_ty() -> Expr {
476    let rn_r = fn_ty(vec_ty(), real_ty());
477    arrow(
478        rn_r.clone(),
479        arrow(list_ty(rn_r), fn_ty(vec_ty(), fn_ty(vec_ty(), real_ty()))),
480    )
481}
482/// `DualFunction : (List Real -> Real) -> (List (List Real -> Real)) -> List Real -> Real`
483/// q(λ) = inf_x L(x, λ) — Lagrange dual function (always concave).
484pub fn dual_function_ty() -> Expr {
485    let rn_r = fn_ty(vec_ty(), real_ty());
486    arrow(
487        rn_r.clone(),
488        arrow(list_ty(rn_r), fn_ty(vec_ty(), real_ty())),
489    )
490}
491/// `WeakDuality : Prop`
492/// Weak duality: d* = sup_λ q(λ) ≤ p* = inf_x f(x).
493pub fn weak_duality_ty() -> Expr {
494    prop()
495}
496/// `StrongDuality : (List Real -> Real) -> (List (List Real -> Real)) -> Prop`
497/// Strong duality: d* = p* (under constraint qualification like Slater's).
498pub fn strong_duality_ty() -> Expr {
499    let rn_r = fn_ty(vec_ty(), real_ty());
500    arrow(rn_r.clone(), arrow(list_ty(rn_r), prop()))
501}
502/// `SlaterCondition : (List (List Real -> Real)) -> (List Real -> Prop) -> Prop`
503/// ∃ x ∈ int(D) with g_i(x) < 0 for all i (strict feasibility).
504pub fn slater_condition_ty() -> Expr {
505    let rn_r = fn_ty(vec_ty(), real_ty());
506    arrow(list_ty(rn_r), arrow(fn_ty(vec_ty(), prop()), prop()))
507}
508/// `KKTConditions : (List Real -> Real) -> (List (List Real -> Real)) -> List Real -> List Real -> Prop`
509/// KKT conditions: stationarity, primal feasibility, dual feasibility, complementary slackness.
510pub fn kkt_conditions_ty() -> Expr {
511    let rn_r = fn_ty(vec_ty(), real_ty());
512    arrow(
513        rn_r.clone(),
514        arrow(list_ty(rn_r), arrow(vec_ty(), arrow(vec_ty(), prop()))),
515    )
516}
517/// `KKTSufficiency : (List Real -> Real) -> (List (List Real -> Real)) -> Prop`
518/// For convex problems, KKT conditions are sufficient for global optimality.
519pub fn kkt_sufficiency_ty() -> Expr {
520    let rn_r = fn_ty(vec_ty(), real_ty());
521    arrow(rn_r.clone(), arrow(list_ty(rn_r), prop()))
522}
523/// `DualityGap : (List Real -> Real) -> (List (List Real -> Real)) -> Real`
524/// Gap = p* - d* ≥ 0.
525pub fn duality_gap_ty() -> Expr {
526    let rn_r = fn_ty(vec_ty(), real_ty());
527    arrow(rn_r.clone(), arrow(list_ty(rn_r), real_ty()))
528}
529/// `FenchelRockafellarDuality : Prop`
530/// inf_x { f(Ax) + g(x) } = - inf_y { f*(y) + g*(-A^T y) } under regularity.
531pub fn fenchel_rockafellar_duality_ty() -> Expr {
532    prop()
533}
534/// `ClarkeSubdifferential : (List Real -> Real) -> List Real -> (List Real -> Prop)`
535/// ∂^C f(x) = conv { lim ∇f(x_k) | x_k → x, x_k ∉ null set }.
536pub fn clarke_subdifferential_ty() -> Expr {
537    arrow(
538        fn_ty(vec_ty(), real_ty()),
539        arrow(vec_ty(), fn_ty(vec_ty(), prop())),
540    )
541}
542/// `ClarkeGeneralisedGradient : (List Real -> Real) -> List Real -> List Real -> Prop`
543/// v ∈ ∂^C f(x): the Clarke generalised directional derivative condition.
544pub fn clarke_generalised_gradient_ty() -> Expr {
545    arrow(
546        fn_ty(vec_ty(), real_ty()),
547        arrow(vec_ty(), arrow(vec_ty(), prop())),
548    )
549}
550/// `IsRegularFunction : (List Real -> Real) -> Prop`
551/// f is Clarke-regular: directional derivative equals Clarke directional derivative.
552pub fn is_regular_function_ty() -> Expr {
553    arrow(fn_ty(vec_ty(), real_ty()), prop())
554}
555/// `NonSmoothChainRule : Prop`
556/// Clarke chain rule: ∂^C(f ∘ g)(x) ⊆ ∂^C f(g(x)) ∘ ∂^C g(x) (in appropriate sense).
557pub fn nonsmooth_chain_rule_ty() -> Expr {
558    prop()
559}
560/// `IsMonotoneOperator : (List Real -> List Real -> Prop) -> Prop`
561/// T is monotone: ⟨Tx - Ty, x - y⟩ ≥ 0 for all x, y.
562pub fn is_monotone_operator_ty() -> Expr {
563    let set_rel = fn_ty(vec_ty(), fn_ty(vec_ty(), prop()));
564    arrow(set_rel, prop())
565}
566/// `IsMaximalMonotone : (List Real -> List Real -> Prop) -> Prop`
567/// T is maximal monotone: no proper monotone extension exists.
568pub fn is_maximal_monotone_ty() -> Expr {
569    let set_rel = fn_ty(vec_ty(), fn_ty(vec_ty(), prop()));
570    arrow(set_rel, prop())
571}
572/// `Resolvent : (List Real -> List Real -> Prop) -> Real -> (List Real -> List Real)`
573/// J_{λT} = (I + λT)^{-1} — resolvent operator of T with parameter λ.
574pub fn resolvent_ty() -> Expr {
575    let set_rel = fn_ty(vec_ty(), fn_ty(vec_ty(), prop()));
576    arrow(set_rel, arrow(real_ty(), fn_ty(vec_ty(), vec_ty())))
577}
578/// `YosidaApproximation : (List Real -> List Real -> Prop) -> Real -> (List Real -> List Real)`
579/// T_λ = (1/λ)(I - J_{λT}) — Yosida approximation (Lipschitz, monotone).
580pub fn yosida_approximation_ty() -> Expr {
581    let set_rel = fn_ty(vec_ty(), fn_ty(vec_ty(), prop()));
582    arrow(set_rel, arrow(real_ty(), fn_ty(vec_ty(), vec_ty())))
583}
584/// `BrouwerFixedPoint : (List Real -> List Real) -> (List Real -> Prop) -> Prop`
585/// Every continuous function f : C → C on a compact convex set has a fixed point.
586pub fn brouwer_fixed_point_ty() -> Expr {
587    arrow(
588        fn_ty(vec_ty(), vec_ty()),
589        arrow(fn_ty(vec_ty(), prop()), prop()),
590    )
591}
592/// `SchauderFixedPoint : (List Real -> List Real) -> (List Real -> Prop) -> Prop`
593/// Infinite-dimensional generalisation: compact continuous map on convex compact set.
594pub fn schauder_fixed_point_ty() -> Expr {
595    arrow(
596        fn_ty(vec_ty(), vec_ty()),
597        arrow(fn_ty(vec_ty(), prop()), prop()),
598    )
599}
600/// `MinimaxTheorem : (List Real -> List Real -> Real) -> (List Real -> Prop) -> (List Real -> Prop) -> Prop`
601/// Rockafellar minimax: min_x max_y f(x,y) = max_y min_x f(x,y) under convex-concave.
602pub fn minimax_theorem_ty() -> Expr {
603    let bilinear = fn_ty(vec_ty(), fn_ty(vec_ty(), real_ty()));
604    let set_ty = fn_ty(vec_ty(), prop());
605    arrow(bilinear, arrow(set_ty.clone(), arrow(set_ty, prop())))
606}
607/// `OperatorSplitting : Prop`
608/// Douglas-Rachford / ADMM splitting convergence for sum of maximal monotone operators.
609pub fn operator_splitting_ty() -> Expr {
610    prop()
611}
612/// `ProjectionOperator : (List Real -> Prop) -> List Real -> List Real`
613/// proj_C(x) = argmin_{y ∈ C} ‖y - x‖.
614pub fn projection_operator_ty() -> Expr {
615    arrow(fn_ty(vec_ty(), prop()), fn_ty(vec_ty(), vec_ty()))
616}
617/// `IsNonexpansive : (List Real -> List Real) -> Prop`
618/// ‖Tx - Ty‖ ≤ ‖x - y‖ for all x, y.
619pub fn is_nonexpansive_ty() -> Expr {
620    arrow(fn_ty(vec_ty(), vec_ty()), prop())
621}
622/// `IsFirmlyNonexpansive : (List Real -> List Real) -> Prop`
623/// ‖Tx - Ty‖² ≤ ⟨Tx - Ty, x - y⟩ for all x, y.
624pub fn is_firmly_nonexpansive_ty() -> Expr {
625    arrow(fn_ty(vec_ty(), vec_ty()), prop())
626}
627/// `AlternatingProjectionConvergence : (List Real -> Prop) -> (List Real -> Prop) -> Prop`
628/// von Neumann / Bregman: alternating projections converge to a point in C ∩ D.
629pub fn alternating_projection_convergence_ty() -> Expr {
630    let set_ty = fn_ty(vec_ty(), prop());
631    arrow(set_ty.clone(), arrow(set_ty, prop()))
632}
633/// Build an [`Environment`] containing convex analysis axioms and theorems.
634pub fn build_convex_analysis_env() -> Environment {
635    let mut env = Environment::new();
636    let axioms: &[(&str, Expr)] = &[
637        ("IsConvexSet", is_convex_set_ty()),
638        ("IsConvexFunction", is_convex_function_ty()),
639        ("IsStrictlyConvexFunction", is_strictly_convex_function_ty()),
640        ("IsStronglyConvexFunction", is_strongly_convex_function_ty()),
641        ("Epigraph", epigraph_ty()),
642        ("LevelSet", level_set_ty()),
643        ("ClosureFunction", closure_function_ty()),
644        ("IsLowerSemicontinuous", is_lsc_ty()),
645        ("IsProperConvex", is_proper_convex_ty()),
646        ("SupportingHyperplane", supporting_hyperplane_ty()),
647        ("SeparatingHyperplane", separating_hyperplane_ty()),
648        ("StrictSeparation", strict_separation_ty()),
649        ("SupportFunction", support_function_ty()),
650        ("GaugeFunction", gauge_function_ty()),
651        ("FenchelConjugate", fenchel_conjugate_ty()),
652        ("Biconjugate", biconjugate_ty()),
653        ("FenchelYoungInequality", fenchel_young_inequality_ty()),
654        ("FenchelYoungEquality", fenchel_young_equality_ty()),
655        ("LegendreFenchelTransform", legendre_fenchel_transform_ty()),
656        ("ConjugateOfSum", conjugate_of_sum_ty()),
657        ("MoreauIdentity", moreau_identity_ty()),
658        ("Subgradient", subgradient_ty()),
659        ("Subdifferential", subdifferential_ty()),
660        ("NormalCone", normal_cone_ty()),
661        ("TangentCone", tangent_cone_ty()),
662        (
663            "OptimalityConditionConvex",
664            optimality_condition_convex_ty(),
665        ),
666        ("SubdiffOfSum", subdiff_of_sum_ty()),
667        ("SubdiffOfComposition", subdiff_of_composition_ty()),
668        ("RecessionCone", recession_cone_ty()),
669        ("RecessionFunction", recession_function_ty()),
670        ("IsCoercive", is_coercive_ty()),
671        ("IsLevelBounded", is_level_bounded_ty()),
672        ("ExistenceOfMinimum", existence_of_minimum_ty()),
673        ("ProximalOperator", proximal_operator_ty()),
674        ("MoreauEnvelope", moreau_envelope_ty()),
675        ("InfConvolution", inf_convolution_ty()),
676        ("ProxFirmlyNonexpansive", prox_firmly_nonexpansive_ty()),
677        ("MoreauEnvelopeSmooth", moreau_envelope_smooth_ty()),
678        (
679            "MoreauEnvelopeConvergence",
680            moreau_envelope_convergence_ty(),
681        ),
682        ("BregmanDivergence", bregman_divergence_ty()),
683        ("BregmanProjection", bregman_projection_ty()),
684        (
685            "BregmanThreePointIdentity",
686            bregman_three_point_identity_ty(),
687        ),
688        ("MirrorDescentStep", mirror_descent_step_ty()),
689        ("ExtremePoint", extreme_point_ty()),
690        ("CaratheodoryTheorem", caratheodory_theorem_ty()),
691        ("KreinMilmanTheorem", krein_milman_theorem_ty()),
692        ("ConvexHull", convex_hull_ty()),
693        ("ClosedConvexHull", closed_convex_hull_ty()),
694        ("FarkasLemma", farkas_lemma_ty()),
695        ("FarkasLemmaMatrix", farkas_lemma_matrix_ty()),
696        ("GordonAlternative", gordon_alternative_ty()),
697        ("TuckerAlternative", tucker_alternative_ty()),
698        ("StiemkeAlternative", stiemke_alternative_ty()),
699        ("IsLipschitzGradient", is_lipschitz_gradient_ty()),
700        (
701            "StrongConvexityQuadraticBound",
702            strong_convexity_quadratic_bound_ty(),
703        ),
704        ("IsSmooth", is_smooth_ty()),
705        ("DescentLemma", descent_lemma_ty()),
706        ("ConvergenceRateGradient", convergence_rate_gradient_ty()),
707        (
708            "LinearConvergenceStrongConvex",
709            linear_convergence_strong_convex_ty(),
710        ),
711        ("IsBarrierFunction", is_barrier_function_ty()),
712        ("IsSelfConcordant", is_self_concordant_ty()),
713        ("LogBarrier", log_barrier_ty()),
714        ("CentralPath", central_path_ty()),
715        ("InteriorPointConvergence", interior_point_convergence_ty()),
716        ("IsSOCConstraint", is_soc_constraint_ty()),
717        ("LorentzCone", lorentz_cone_ty()),
718        ("PositiveSemidefiniteCone", psd_cone_ty()),
719        ("SOCPDuality", socp_duality_ty()),
720        ("SDPDuality", sdp_duality_ty()),
721        ("LagrangianFunction", lagrangian_function_ty()),
722        ("DualFunction", dual_function_ty()),
723        ("WeakDuality", weak_duality_ty()),
724        ("StrongDuality", strong_duality_ty()),
725        ("SlaterCondition", slater_condition_ty()),
726        ("KKTConditions", kkt_conditions_ty()),
727        ("KKTSufficiency", kkt_sufficiency_ty()),
728        ("DualityGap", duality_gap_ty()),
729        (
730            "FenchelRockafellarDuality",
731            fenchel_rockafellar_duality_ty(),
732        ),
733        ("ClarkeSubdifferential", clarke_subdifferential_ty()),
734        (
735            "ClarkeGeneralisedGradient",
736            clarke_generalised_gradient_ty(),
737        ),
738        ("IsRegularFunction", is_regular_function_ty()),
739        ("NonSmoothChainRule", nonsmooth_chain_rule_ty()),
740        ("IsMonotoneOperator", is_monotone_operator_ty()),
741        ("IsMaximalMonotone", is_maximal_monotone_ty()),
742        ("Resolvent", resolvent_ty()),
743        ("YosidaApproximation", yosida_approximation_ty()),
744        ("BrouwerFixedPoint", brouwer_fixed_point_ty()),
745        ("SchauderFixedPoint", schauder_fixed_point_ty()),
746        ("MinimaxTheorem", minimax_theorem_ty()),
747        ("OperatorSplitting", operator_splitting_ty()),
748        ("ProjectionOperator", projection_operator_ty()),
749        ("IsNonexpansive", is_nonexpansive_ty()),
750        ("IsFirmlyNonexpansive", is_firmly_nonexpansive_ty()),
751        (
752            "AlternatingProjectionConvergence",
753            alternating_projection_convergence_ty(),
754        ),
755    ];
756    for (name, ty) in axioms {
757        let decl = Declaration::Axiom {
758            name: Name::str(*name),
759            univ_params: vec![],
760            ty: ty.clone(),
761        };
762        let _ = env.add(decl);
763    }
764    env
765}
766/// Computes a subgradient of a convex function at a given point via finite differences.
767pub fn compute_subgradient(f: &ConvexFunction, x: &[f64]) -> Vec<f64> {
768    f.gradient(x)
769}
770/// Checks the subgradient inequality: f(y) ≥ f(x) + ⟨g, y-x⟩.
771pub fn check_subgradient_inequality(f: &ConvexFunction, x: &[f64], g: &[f64], y: &[f64]) -> bool {
772    let fx = f.eval(x);
773    let fy = f.eval(y);
774    let inner: f64 = g
775        .iter()
776        .zip(y.iter().zip(x.iter()))
777        .map(|(gi, (yi, xi))| gi * (yi - xi))
778        .sum();
779    fy >= fx + inner - 1e-9
780}
781/// Checks v ∈ N_C(x): ⟨v, y - x⟩ ≤ ε for all sampled y in C.
782pub fn check_normal_cone_membership(v: &[f64], x: &[f64], set_points: &[Vec<f64>]) -> bool {
783    for y in set_points {
784        let dot: f64 = v
785            .iter()
786            .zip(y.iter().zip(x.iter()))
787            .map(|(vi, (yi, xi))| vi * (yi - xi))
788            .sum();
789        if dot > 1e-9 {
790            return false;
791        }
792    }
793    true
794}
795/// Compute prox_{λf}(v) = argmin_x { f(x) + 1/(2λ) ‖x-v‖² } using gradient descent.
796pub fn proximal_operator(f: &ConvexFunction, v: &[f64], cfg: &ProxConfig) -> Vec<f64> {
797    let n = v.len();
798    let mut x = v.to_vec();
799    let inv_lambda = 1.0 / cfg.lambda;
800    for _ in 0..cfg.max_iter {
801        let grad_f = f.gradient(&x);
802        let mut new_x = vec![0.0; n];
803        let mut delta_sq = 0.0;
804        for i in 0..n {
805            let g = grad_f[i] + inv_lambda * (x[i] - v[i]);
806            let dx = cfg.step_size * g;
807            new_x[i] = x[i] - dx;
808            delta_sq += dx * dx;
809        }
810        x = new_x;
811        if delta_sq.sqrt() < cfg.tol {
812            break;
813        }
814    }
815    x
816}
817/// Compute Moreau envelope e_{λf}(x) = inf_y { f(y) + 1/(2λ)‖x-y‖² }.
818pub fn moreau_envelope(f: &ConvexFunction, x: &[f64], cfg: &ProxConfig) -> f64 {
819    let p = proximal_operator(f, x, cfg);
820    let dist_sq: f64 = x
821        .iter()
822        .zip(p.iter())
823        .map(|(xi, pi)| (xi - pi).powi(2))
824        .sum();
825    f.eval(&p) + dist_sq / (2.0 * cfg.lambda)
826}
827/// Compute inf-convolution (f □ g)(x) = inf_y { f(y) + g(x-y) } by grid search in 1D.
828pub fn inf_convolution_1d(
829    f: fn(f64) -> f64,
830    g: fn(f64) -> f64,
831    x: f64,
832    steps: usize,
833    radius: f64,
834) -> f64 {
835    let step = 2.0 * radius / steps as f64;
836    let mut best = f64::INFINITY;
837    for k in 0..=steps {
838        let y = -radius + k as f64 * step;
839        let val = f(y) + g(x - y);
840        if val < best {
841            best = val;
842        }
843    }
844    best
845}
846/// Bregman divergence D_f(x, y) = f(x) - f(y) - ⟨∇f(y), x-y⟩.
847pub fn bregman_divergence(f: &ConvexFunction, x: &[f64], y: &[f64]) -> f64 {
848    let fy = f.eval(y);
849    let fx = f.eval(x);
850    let grad_y = f.gradient(y);
851    let inner: f64 = grad_y
852        .iter()
853        .zip(x.iter().zip(y.iter()))
854        .map(|(g, (xi, yi))| g * (xi - yi))
855        .sum();
856    fx - fy - inner
857}
858/// Check three-point identity (approximately):
859/// D_f(x, z) ≈ D_f(x, y) + D_f(y, z) + ⟨∇f(y) - ∇f(z), x - y⟩.
860pub fn check_three_point_identity(f: &ConvexFunction, x: &[f64], y: &[f64], z: &[f64]) -> bool {
861    let lhs = bregman_divergence(f, x, z);
862    let d_xy = bregman_divergence(f, x, y);
863    let d_yz = bregman_divergence(f, y, z);
864    let grad_y = f.gradient(y);
865    let grad_z = f.gradient(z);
866    let inner: f64 = grad_y
867        .iter()
868        .zip(grad_z.iter())
869        .zip(x.iter().zip(y.iter()))
870        .map(|((gy, gz), (xi, yi))| (gy - gz) * (xi - yi))
871        .sum();
872    let rhs = d_xy + d_yz + inner;
873    (lhs - rhs).abs() < 1e-5
874}
875/// Compute a separating hyperplane between two finite point sets using the midpoint normal.
876pub fn compute_separating_hyperplane(
877    a_points: &[Vec<f64>],
878    b_points: &[Vec<f64>],
879) -> Option<Hyperplane> {
880    if a_points.is_empty() || b_points.is_empty() {
881        return None;
882    }
883    let n = a_points[0].len();
884    let ca: Vec<f64> = (0..n)
885        .map(|i| a_points.iter().map(|p| p[i]).sum::<f64>() / a_points.len() as f64)
886        .collect();
887    let cb: Vec<f64> = (0..n)
888        .map(|i| b_points.iter().map(|p| p[i]).sum::<f64>() / b_points.len() as f64)
889        .collect();
890    let normal: Vec<f64> = ca.iter().zip(cb.iter()).map(|(a, b)| a - b).collect();
891    let norm = normal.iter().map(|x| x * x).sum::<f64>().sqrt();
892    if norm < 1e-12 {
893        return None;
894    }
895    let unit_normal: Vec<f64> = normal.iter().map(|x| x / norm).collect();
896    let mid: Vec<f64> = ca
897        .iter()
898        .zip(cb.iter())
899        .map(|(a, b)| (a + b) / 2.0)
900        .collect();
901    let offset: f64 = unit_normal.iter().zip(mid.iter()).map(|(a, b)| a * b).sum();
902    Some(Hyperplane::new(unit_normal, offset))
903}
904/// Compute σ_C(y) = max_{x ∈ C} ⟨y, x⟩ over a finite set of points.
905pub fn support_function(c_points: &[Vec<f64>], y: &[f64]) -> f64 {
906    c_points
907        .iter()
908        .map(|x| x.iter().zip(y.iter()).map(|(xi, yi)| xi * yi).sum::<f64>())
909        .fold(f64::NEG_INFINITY, f64::max)
910}
911/// Compute gauge function γ_C(x) = inf { t ≥ 0 | x ∈ t·C } via binary search.
912pub fn gauge_function(c_points: &[Vec<f64>], x: &[f64]) -> f64 {
913    let norm_x: f64 = x.iter().map(|xi| xi * xi).sum::<f64>().sqrt();
914    if norm_x < 1e-15 {
915        return 0.0;
916    }
917    let mut lo = 0.0_f64;
918    let mut hi = norm_x * 100.0 + 1.0;
919    for _ in 0..60 {
920        let mid = (lo + hi) / 2.0;
921        let scaled: Vec<f64> = x.iter().map(|xi| xi / mid.max(1e-15)).collect();
922        let s = support_function(c_points, &scaled);
923        let s_x = support_function(c_points, x);
924        if s_x <= mid * s + 1e-9 {
925            hi = mid;
926        } else {
927            lo = mid;
928        }
929    }
930    hi
931}
932#[cfg(test)]
933mod tests {
934    use super::*;
935    fn sq_func(x: &[f64]) -> f64 {
936        x.iter().map(|xi| xi * xi).sum::<f64>() * 0.5
937    }
938    fn sq_grad(x: &[f64]) -> Vec<f64> {
939        x.to_vec()
940    }
941    #[test]
942    fn test_build_env_keys() {
943        let env = build_convex_analysis_env();
944        assert!(env.get(&Name::str("IsConvexFunction")).is_some());
945        assert!(env.get(&Name::str("FenchelConjugate")).is_some());
946        assert!(env.get(&Name::str("Subdifferential")).is_some());
947        assert!(env.get(&Name::str("NormalCone")).is_some());
948        assert!(env.get(&Name::str("ProximalOperator")).is_some());
949        assert!(env.get(&Name::str("MoreauEnvelope")).is_some());
950        assert!(env.get(&Name::str("BregmanDivergence")).is_some());
951        assert!(env.get(&Name::str("InfConvolution")).is_some());
952    }
953    #[test]
954    fn test_epigraph_membership() {
955        let f = ConvexFunction::new("sq", sq_func, Some(sq_grad));
956        assert!(f.in_epigraph(&[1.0, 1.0], 1.5));
957        assert!(!f.in_epigraph(&[1.0, 1.0], 0.5));
958    }
959    #[test]
960    fn test_level_set_membership() {
961        let f = ConvexFunction::new("sq", sq_func, Some(sq_grad));
962        assert!(f.in_level_set(&[0.5, 0.5], 1.0));
963        assert!(!f.in_level_set(&[2.0, 0.0], 1.0));
964    }
965    #[test]
966    fn test_subgradient_inequality() {
967        let f = ConvexFunction::new("sq", sq_func, Some(sq_grad));
968        let x = vec![1.0];
969        let g = f.gradient(&x);
970        let y = vec![3.0];
971        assert!(check_subgradient_inequality(&f, &x, &g, &y));
972        let y2 = vec![-2.0];
973        assert!(check_subgradient_inequality(&f, &x, &g, &y2));
974    }
975    #[test]
976    fn test_proximal_operator_quadratic() {
977        let f = ConvexFunction::new("sq", sq_func, Some(sq_grad));
978        let cfg = ProxConfig::new(1.0);
979        let v = vec![2.0];
980        let p = proximal_operator(&f, &v, &cfg);
981        assert!((p[0] - 1.0).abs() < 1e-3, "prox = {}", p[0]);
982    }
983    #[test]
984    fn test_moreau_envelope_quadratic() {
985        let f = ConvexFunction::new("sq", sq_func, Some(sq_grad));
986        let cfg = ProxConfig::new(1.0);
987        let x = vec![2.0];
988        let env_val = moreau_envelope(&f, &x, &cfg);
989        let expected = 0.5 * 4.0 / 2.0;
990        assert!((env_val - expected).abs() < 1e-3, "env = {}", env_val);
991    }
992    #[test]
993    fn test_bregman_divergence_quadratic() {
994        let f = ConvexFunction::new("sq", sq_func, Some(sq_grad));
995        let x = vec![3.0, 0.0];
996        let y = vec![1.0, 0.0];
997        let d = bregman_divergence(&f, &x, &y);
998        assert!((d - 2.0).abs() < 1e-6, "D_f = {d}");
999    }
1000    #[test]
1001    fn test_bregman_three_point_identity() {
1002        let f = ConvexFunction::new("sq", sq_func, Some(sq_grad));
1003        let x = vec![3.0];
1004        let y = vec![1.0];
1005        let z = vec![-1.0];
1006        assert!(check_three_point_identity(&f, &x, &y, &z));
1007    }
1008    #[test]
1009    fn test_separating_hyperplane() {
1010        let a = vec![vec![0.0, 0.0], vec![0.5, 0.0]];
1011        let b = vec![vec![2.0, 0.0], vec![2.5, 0.0]];
1012        let hp = compute_separating_hyperplane(&a, &b).expect("operation should succeed");
1013        assert!(hp.separates(&a, &b), "hyperplane should separate A from B");
1014    }
1015    #[test]
1016    fn test_support_function_unit_square() {
1017        let c = vec![
1018            vec![0.0, 0.0],
1019            vec![1.0, 0.0],
1020            vec![1.0, 1.0],
1021            vec![0.0, 1.0],
1022        ];
1023        let s1 = support_function(&c, &[1.0, 0.0]);
1024        assert!((s1 - 1.0).abs() < 1e-12, "σ(1,0) = 1");
1025        let s2 = support_function(&c, &[1.0, 1.0]);
1026        assert!((s2 - 2.0).abs() < 1e-12, "σ(1,1) = 2");
1027    }
1028    #[test]
1029    fn test_inf_convolution_1d() {
1030        let f = |x: f64| 0.5 * x * x;
1031        let g = |x: f64| 0.5 * x * x;
1032        let x = 4.0;
1033        let val = inf_convolution_1d(f, g, x, 1000, 10.0);
1034        let expected = 0.25 * x * x;
1035        assert!(
1036            (val - expected).abs() < 0.05,
1037            "inf-conv = {val}, expected {expected}"
1038        );
1039    }
1040    #[test]
1041    fn test_build_env_new_axioms() {
1042        let env = build_convex_analysis_env();
1043        assert!(env.get(&Name::str("ExtremePoint")).is_some());
1044        assert!(env.get(&Name::str("CaratheodoryTheorem")).is_some());
1045        assert!(env.get(&Name::str("KreinMilmanTheorem")).is_some());
1046        assert!(env.get(&Name::str("ConvexHull")).is_some());
1047        assert!(env.get(&Name::str("FarkasLemma")).is_some());
1048        assert!(env.get(&Name::str("GordonAlternative")).is_some());
1049        assert!(env.get(&Name::str("IsLipschitzGradient")).is_some());
1050        assert!(env.get(&Name::str("IsSmooth")).is_some());
1051        assert!(env.get(&Name::str("DescentLemma")).is_some());
1052        assert!(env.get(&Name::str("IsBarrierFunction")).is_some());
1053        assert!(env.get(&Name::str("IsSelfConcordant")).is_some());
1054        assert!(env.get(&Name::str("CentralPath")).is_some());
1055        assert!(env.get(&Name::str("LorentzCone")).is_some());
1056        assert!(env.get(&Name::str("PositiveSemidefiniteCone")).is_some());
1057        assert!(env.get(&Name::str("SOCPDuality")).is_some());
1058        assert!(env.get(&Name::str("KKTConditions")).is_some());
1059        assert!(env.get(&Name::str("SlaterCondition")).is_some());
1060        assert!(env.get(&Name::str("StrongDuality")).is_some());
1061        assert!(env.get(&Name::str("FenchelRockafellarDuality")).is_some());
1062        assert!(env.get(&Name::str("ClarkeSubdifferential")).is_some());
1063        assert!(env.get(&Name::str("IsRegularFunction")).is_some());
1064        assert!(env.get(&Name::str("IsMonotoneOperator")).is_some());
1065        assert!(env.get(&Name::str("Resolvent")).is_some());
1066        assert!(env.get(&Name::str("BrouwerFixedPoint")).is_some());
1067        assert!(env.get(&Name::str("MinimaxTheorem")).is_some());
1068        assert!(env.get(&Name::str("ProjectionOperator")).is_some());
1069        assert!(env
1070            .get(&Name::str("AlternatingProjectionConvergence"))
1071            .is_some());
1072    }
1073    #[test]
1074    fn test_subgradient_method_quadratic() {
1075        let f = ConvexFunction::new("sq", sq_func, Some(sq_grad));
1076        let method = SubgradientMethod::new(StepSchedule::DiminishingSqrt(1.0), 500);
1077        let (best, _history) = method.run(&f, &[5.0]);
1078        assert!(
1079            best[0].abs() < 0.5,
1080            "should converge near 0, got {}",
1081            best[0]
1082        );
1083    }
1084    #[test]
1085    fn test_proximal_point_algorithm_quadratic() {
1086        let f = ConvexFunction::new("sq", sq_func, Some(sq_grad));
1087        let ppa = ProximalPointAlgorithm::constant(1.0, 50, 1e-6);
1088        let iterates = ppa.run(&f, &[4.0]);
1089        let last = iterates.last().expect("last should succeed");
1090        assert!(
1091            last[0].abs() < 0.1,
1092            "PPA should converge near 0, got {}",
1093            last[0]
1094        );
1095    }
1096    #[test]
1097    fn test_fenchel_conjugate_squared_norm() {
1098        let ev = FenchelConjugateEvaluator::new(FunctionClass::SquaredNorm);
1099        let y = vec![3.0, 4.0];
1100        let conj = ev.eval(&y);
1101        assert!((conj - 12.5).abs() < 1e-10, "f*(3,4) = 12.5, got {conj}");
1102    }
1103    #[test]
1104    fn test_fenchel_young_inequality() {
1105        let ev = FenchelConjugateEvaluator::new(FunctionClass::SquaredNorm);
1106        let x = vec![1.0, 2.0];
1107        let y = vec![3.0, 0.5];
1108        let fx = 0.5 * (1.0_f64 + 4.0_f64);
1109        assert!(ev.check_fenchel_young(&x, &y, fx));
1110    }
1111    #[test]
1112    fn test_fenchel_conjugate_negative_entropy() {
1113        let ev = FenchelConjugateEvaluator::new(FunctionClass::NegativeEntropy);
1114        let y = vec![1.0];
1115        let val = ev.eval(&y);
1116        assert!((val - 1.0).abs() < 1e-10, "f*(1) = exp(0) = 1, got {val}");
1117    }
1118    #[test]
1119    fn test_fenchel_conjugate_box_indicator() {
1120        let ev = FenchelConjugateEvaluator::new(FunctionClass::BoxIndicator { lo: -1.0, hi: 1.0 });
1121        let y = vec![2.0, -3.0];
1122        let val = ev.eval(&y);
1123        assert!((val - 5.0).abs() < 1e-10, "f*(2,-3) = 5, got {val}");
1124    }
1125    #[test]
1126    fn test_separating_hyperplane_finder() {
1127        let a = vec![vec![0.0, 0.0], vec![0.0, 1.0]];
1128        let b = vec![vec![3.0, 0.0], vec![3.0, 1.0]];
1129        let finder = SeparatingHyperplaneFinder::new(0.01, 200, 1e-4);
1130        let hp = finder.find(&a, &b);
1131        assert!(hp.is_some(), "should find a separating hyperplane");
1132        let hp = hp.expect("hp should be valid");
1133        assert!(hp.separates(&a, &b), "hyperplane should separate A and B");
1134    }
1135    #[test]
1136    fn test_alternating_projections_intersection() {
1137        let proj_a = |x: &[f64]| vec![x[0].max(1.0)];
1138        let proj_b = |x: &[f64]| vec![x[0].min(2.0)];
1139        let solver = AlternatingProjectionSolver::new(100, 1e-8);
1140        let iters = solver.run(proj_a, proj_b, &[5.0]);
1141        let last = iters.last().expect("last should succeed");
1142        assert!(
1143            (last[0] - 2.0).abs() < 1e-6,
1144            "should converge to 2, got {}",
1145            last[0]
1146        );
1147    }
1148}
1149#[cfg(test)]
1150mod convex_ext_tests {
1151    use super::*;
1152    #[test]
1153    fn test_convex_cone() {
1154        let soc = ConvexCone::second_order_cone(3);
1155        assert!(soc.is_pointed);
1156        assert!(soc.is_closed);
1157        assert!(!soc.dual_cone_description().is_empty());
1158    }
1159    #[test]
1160    fn test_fenchel_conjugate() {
1161        let fc = FenchelConjugate::new("||x||^2/2");
1162        assert!(!fc.definition().is_empty());
1163        assert!(!fc.fenchel_inequality().is_empty());
1164    }
1165    #[test]
1166    fn test_proximal_operator() {
1167        let prox = ProximalOperator::new("||x||_1", 0.1);
1168        assert!(!prox.definition().is_empty());
1169        assert!(!prox.for_indicator_function().is_empty());
1170    }
1171    #[test]
1172    fn test_admm() {
1173        let admm = AdmmSolver::new(1.0, 1000, 1e-4);
1174        assert!(!admm.convergence_description().is_empty());
1175    }
1176    #[test]
1177    fn test_convex_program() {
1178        let lp = ConvexProgram::new("diet", 100, 50, ConvexProblemClass::Lp);
1179        assert!(lp.is_lp());
1180        assert!(lp.strong_duality_holds());
1181    }
1182}
1183#[cfg(test)]
1184mod duality_ext_tests {
1185    use super::*;
1186    #[test]
1187    fn test_lagrangian_duality() {
1188        let ld = LagrangianDuality::new("min f(x)", "max g(lambda)", 0.0);
1189        assert!(ld.strong_duality());
1190        assert_eq!(ld.kkt_conditions().len(), 4);
1191    }
1192    #[test]
1193    fn test_mirror_descent() {
1194        let md = MirrorDescent::with_entropy(0.01, 1000);
1195        let rate = md.convergence_rate(1.0, 1.0);
1196        assert!(rate > 0.0 && rate < 1.0);
1197    }
1198}
1199#[cfg(test)]
1200mod epi_level_tests {
1201    use super::*;
1202    #[test]
1203    fn test_epigraph() {
1204        let epi = Epigraph::new("||x||^2");
1205        assert!(epi.f_convex_iff_epi_convex());
1206        assert!(!epi.definition().is_empty());
1207    }
1208    #[test]
1209    fn test_sublevel_set() {
1210        let sl = SublevelSet::new("exp(x)", 2.0);
1211        assert!(sl.is_convex_if_f_quasiconvex());
1212    }
1213    #[test]
1214    fn test_recession_cone() {
1215        let rc = RecessionCone::new("R^n+");
1216        assert!(rc.compact_iff_trivial_recession());
1217    }
1218}
1219#[cfg(test)]
1220mod tests_convex_analysis_ext {
1221    use super::*;
1222    #[test]
1223    fn test_convex_conjugate() {
1224        let quad = ConvexConjugate::quadratic_conjugate();
1225        assert!(quad.biconjugate_equals_f());
1226        let fm = quad.fenchel_moreau_theorem();
1227        assert!(fm.contains("Fenchel-Moreau"));
1228        let yf = quad.young_fenchel_inequality();
1229        assert!(yf.contains("Young-Fenchel"));
1230        let ind = ConvexConjugate::indicator_conjugate("C");
1231        assert!(ind.conjugate_name.contains("h_"));
1232    }
1233    #[test]
1234    fn test_fenchel_duality() {
1235        let strong = FenchelDualityPair::new("min f", "max -f*", 0.0);
1236        assert!(strong.strong_duality_holds);
1237        assert!(strong.slater_condition_met());
1238        let desc = strong.duality_gap_description();
1239        assert!(desc.contains("Strong duality"));
1240    }
1241    #[test]
1242    fn test_proximal_l1() {
1243        let prox = ProximalOperatorNew::l1_norm(0.5);
1244        assert!(prox.has_closed_form);
1245        let formula = prox.proximal_point_formula();
1246        assert!(formula.contains("prox_"));
1247        let moreau = prox.moreau_decomposition();
1248        assert!(moreau.contains("Moreau"));
1249    }
1250    #[test]
1251    fn test_admm() {
1252        let mut admm = ADMMData::new(1.0);
1253        assert!(!admm.has_converged(1e-6));
1254        admm.update_residuals(1e-8, 1e-8);
1255        assert!(admm.has_converged(1e-6));
1256        let desc = admm.admm_update_description();
1257        assert!(desc.contains("ADMM"));
1258        let conv = admm.convergence_guarantee();
1259        assert!(conv.contains("ADMM"));
1260    }
1261    #[test]
1262    fn test_variational_inequality() {
1263        let vi = VariationalInequality::new("F", "C", true);
1264        let form = vi.vi_formulation();
1265        assert!(form.contains("F(x*)"));
1266        let stamp = vi.stampacchia_existence();
1267        assert!(stamp.contains("Stampacchia"));
1268        let svi = VariationalInequality::strongly_monotone("G", "K", 0.5);
1269        assert!(svi.unique_solution_exists());
1270    }
1271    #[test]
1272    fn test_extragradient() {
1273        let mut eg = ExtragradientMethod::new(0.01, "F");
1274        assert!(eg.convergence_condition(50.0));
1275        assert!(!eg.convergence_condition(200.0));
1276        eg.do_step(0.5);
1277        assert_eq!(eg.iterations, 1);
1278        let desc = eg.korpelevich_step_description();
1279        assert!(desc.contains("Korpelevich"));
1280    }
1281}