module public import Causalean.PO.ID.Partial.LP.ConicDuality public import Mathlib.Algebra.BigOperators.Group.Finset.Basic public import Mathlib.Algebra.BigOperators.Group.Finset.Defs public import Mathlib.Algebra.BigOperators.Group.Finset.Piecewise public import Mathlib.Algebra.BigOperators.Group.Finset.Sigma public import Mathlib.Algebra.BigOperators.Ring.Finset public import Mathlib.Algebra.Field.Defs public import Mathlib.Algebra.Group.Action.Defs public import Mathlib.Algebra.Group.Basic public import Mathlib.Algebra.Group.Defs public import Mathlib.Algebra.Group.Hom.Defs public import Mathlib.Algebra.Group.Pi.Basic public import Mathlib.Algebra.Group.Submonoid.Defs public import Mathlib.Algebra.Group.Subsemigroup.Defs public import Mathlib.Algebra.GroupWithZero.Action.Defs public import Mathlib.Algebra.GroupWithZero.Defs public import Mathlib.Algebra.Module.Defs public import Mathlib.Algebra.Module.Equiv.Defs public import Mathlib.Algebra.Module.LinearMap.Defs public import Mathlib.Algebra.Module.Pi public import Mathlib.Algebra.Module.Submodule.Defs public import Mathlib.Algebra.Module.Submodule.Equiv public import Mathlib.Algebra.Notation.Defs public import Mathlib.Algebra.Notation.Pi.Defs public import Mathlib.Algebra.Order.Archimedean.Basic public import Mathlib.Algebra.Order.Archimedean.Real.Basic public import Mathlib.Algebra.Order.BigOperators.Group.Finset public import Mathlib.Algebra.Order.Group.Unbundled.Abs public import Mathlib.Algebra.Order.Group.Unbundled.Basic public import Mathlib.Algebra.Order.GroupWithZero.Basic public import Mathlib.Algebra.Order.Monoid.Defs public import Mathlib.Algebra.Order.Monoid.Unbundled.Basic public import Mathlib.Algebra.Order.Monoid.Unbundled.Defs public import Mathlib.Algebra.Order.Nonneg.Basic public import Mathlib.Algebra.Order.Nonneg.Module public import Mathlib.Algebra.Order.Ring.Defs public import Mathlib.Algebra.Order.Ring.Nat public import Mathlib.Algebra.Ring.CompTypeclasses public import Mathlib.Algebra.Ring.Defs public import Mathlib.Algebra.Ring.Hom.Defs public import Mathlib.Algebra.Ring.Nat public import Mathlib.Algebra.Star.Basic public import Mathlib.Analysis.Convex.Cone.Basic public import Mathlib.Analysis.Convex.Cone.InnerDual public import Mathlib.Analysis.InnerProductSpace.Adjoint public import Mathlib.Analysis.InnerProductSpace.Basic public import Mathlib.Analysis.InnerProductSpace.Continuous public import Mathlib.Analysis.InnerProductSpace.Defs public import Mathlib.Analysis.InnerProductSpace.PiL2 public import Mathlib.Analysis.Normed.Field.Basic public import Mathlib.Analysis.Normed.Group.Defs public import Mathlib.Analysis.Normed.Group.Real public import Mathlib.Analysis.Normed.Group.Uniform public import Mathlib.Analysis.Normed.Lp.PiLp public import Mathlib.Analysis.Normed.Lp.WithLp public import Mathlib.Analysis.Normed.Module.Basic public import Mathlib.Analysis.Normed.Operator.Basic public import Mathlib.Analysis.Normed.Operator.LinearIsometry public import Mathlib.Analysis.Normed.Order.Lattice public import Mathlib.Analysis.Normed.Ring.Basic public import Mathlib.Analysis.RCLike.Basic public import Mathlib.Analysis.RCLike.Lemmas public import Mathlib.Data.Countable.Basic public import Mathlib.Data.Countable.Defs public import Mathlib.Data.ENNReal.Basic public import Mathlib.Data.Finite.Prod public import Mathlib.Data.Finset.Defs public import Mathlib.Data.Fintype.BigOperators public import Mathlib.Data.Fintype.Defs public import Mathlib.Data.Fintype.EquivFin public import Mathlib.Data.Fintype.Prod public import Mathlib.Data.Fintype.Sum public import Mathlib.Data.FunLike.Basic public import Mathlib.Data.FunLike.Equiv public import Mathlib.Data.Int.Cast.Defs public import Mathlib.Data.Matrix.Mul public import Mathlib.Data.NNReal.Defs public import Mathlib.Data.Nat.Basic public import Mathlib.Data.Nat.Cast.Defs public import Mathlib.Data.Nat.Init public import Mathlib.Data.Real.Basic public import Mathlib.Data.Real.Star public import Mathlib.Data.Set.Basic public import Mathlib.Data.Set.Defs public import Mathlib.Data.Set.Image public import Mathlib.Data.Set.Insert public import Mathlib.Data.Set.Lattice public import Mathlib.Data.Set.Operations public import Mathlib.Data.SetLike.Basic public import Mathlib.LinearAlgebra.DFinsupp public import Mathlib.LinearAlgebra.FiniteDimensional.Defs public import Mathlib.Logic.Basic public import Mathlib.Logic.Equiv.Defs public import Mathlib.Order.Basic public import Mathlib.Order.Bounds.Defs public import Mathlib.Order.ConditionallyCompleteLattice.Basic public import Mathlib.Order.ConditionallyCompleteLattice.Defs public import Mathlib.Order.ConditionallyCompletePartialOrder.Defs public import Mathlib.Order.Defs.LinearOrder public import Mathlib.Order.Defs.PartialOrder public import Mathlib.Order.Filter.AtTopBot.Basic public import Mathlib.Order.Filter.AtTopBot.Defs public import Mathlib.Order.Filter.AtTopBot.Tendsto public import Mathlib.Order.Filter.Defs public import Mathlib.Order.Filter.Extr public import Mathlib.Order.Filter.Tendsto public import Mathlib.Order.Interval.Set.Defs public import Mathlib.Order.Lattice public import Mathlib.Order.Monotone.Defs public import Mathlib.Order.RelClasses public import Mathlib.Order.SetNotation public import Mathlib.Tactic.Linarith.Lemmas public import Mathlib.Tactic.NormNum.Basic public import Mathlib.Tactic.NormNum.Ineq public import Mathlib.Tactic.NormNum.Result public import Mathlib.Tactic.Ring.Basic public import Mathlib.Tactic.Ring.Common public import Mathlib.Topology.Algebra.Group.Defs public import Mathlib.Topology.Algebra.Module.ClosedSubmodule public import Mathlib.Topology.Algebra.Module.ContinuousLinearMap.Basic public import Mathlib.Topology.Algebra.Module.Equiv public import Mathlib.Topology.Algebra.Module.FiniteDimension public import Mathlib.Topology.Algebra.Monoid public import Mathlib.Topology.Algebra.Monoid.Defs public import Mathlib.Topology.Algebra.MulAction public import Mathlib.Topology.Algebra.Ring.Basic public import Mathlib.Topology.Algebra.Ring.Real public import Mathlib.Topology.Algebra.UniformMulAction public import Mathlib.Topology.Bases public import Mathlib.Topology.Basic public import Mathlib.Topology.Compactness.Compact public import Mathlib.Topology.Constructions public import Mathlib.Topology.Constructions.SumProd public import Mathlib.Topology.Continuous public import Mathlib.Topology.ContinuousOn public import Mathlib.Topology.Defs.Basic public import Mathlib.Topology.Defs.Filter public import Mathlib.Topology.Defs.Sequences public import Mathlib.Topology.EMetricSpace.Defs public import Mathlib.Topology.MetricSpace.Algebra public import Mathlib.Topology.MetricSpace.Basic public import Mathlib.Topology.MetricSpace.Pseudo.Defs public import Mathlib.Topology.MetricSpace.Pseudo.Lemmas public import Mathlib.Topology.Metrizable.Basic public import Mathlib.Topology.Metrizable.Uniformity public import Mathlib.Topology.Order.Compact public import Mathlib.Topology.Order.OrderClosed public import Mathlib.Topology.Separation.Hausdorff public import Mathlib.Topology.Sequences public import Mathlib.Topology.UniformSpace.Defs public import Mathlib.Topology.UniformSpace.Real /-! # Finite transport linear programs The transportation problem between two finite mass vectors with a general cost matrix: feasible plans, existence of a cost-minimizing plan, and an approximate dual certificate (row and column potentials satisfying every dual constraint whose objective is within any tolerance of the optimal cost). Slot weights may vanish and slot types are not assumed inhabited. The dual certificate is obtained from Causalean's conic strong duality (`Causalean.PartialID.ConicProgram.strong_duality_zero_gap`); the closedness of the augmented image comes from row sums bounding nonnegative plans. Used by the finite-support step of Kantorovich–Rubinstein duality (`Toolkit.Foundations.MeasureTheory.KantorovichRubinstein`). -/ @[expose] public section namespace Toolkit.Foundations.Optimization.Convex open scoped BigOperators NNReal RealInnerProductSpace open Filter Topology variable {ι κ : Type*} [Fintype ι] [Fintype κ] /-- Given [nonnegative row masses](hyp:a) and [nonnegative column masses](hyp:b), a [matrix](hyp:p) indexed by [finite row and column index sets](hyp:ι,κ) is a [finite transport plan](goal) when [all its entries are nonnegative](step:1), [its row sums equal the row masses](step:2), and [its column sums equal the column masses](step:3). -/ def IsFiniteTransportPlan (a : ι → ℝ≥0) (b : κ → ℝ≥0) (p : ι → κ → ℝ) : Prop := (∀ i j, 0 ≤ p i j) ∧ (∀ i, ∑ j, p i j = (a i : ℝ)) ∧ (∀ j, ∑ i, p i j = (b j : ℝ)) /-- For [row masses summing to one](hyp:ha) and [column masses summing to one](hyp:hb), and [any finite cost matrix](hyp:c), [some transport plan between the two mass vectors has total cost no larger than that of any other transport plan](goal). -/ theorem exists_minimal_finiteTransportPlan (a : ι → ℝ≥0) (b : κ → ℝ≥0) (ha : ∑ i, (a i : ℝ) = 1) (hb : ∑ j, (b j : ℝ) = 1) (c : ι → κ → ℝ) : ∃ p, IsFiniteTransportPlan a b p ∧ ∀ q, IsFiniteTransportPlan a b q → (∑ i, ∑ j, c i j * p i j) ≤ ∑ i, ∑ j, c i j * q i j := by -- Product weights give feasibility. Row masses bound every feasible entry. -- The feasible set is closed and bounded in finite-dimensional Euclidean -- coordinates, so the continuous objective attains its minimum. classical let S : Set (ι → κ → ℝ) := {p | IsFiniteTransportPlan a b p} have hne : S.Nonempty := by refine ⟨fun i j => (a i : ℝ) * (b j : ℝ), ?_, ?_, ?_⟩ · intro i j positivity · intro i simp only [← Finset.mul_sum, hb, mul_one] · intro j simp only [← Finset.sum_mul, ha, one_mul] have hclosed : IsClosed S := by have heq : S = (⋂ i, ⋂ j, {p | 0 ≤ p i j}) ∩ ((⋂ i, {p | ∑ j, p i j = (a i : ℝ)}) ∩ (⋂ j, {p | ∑ i, p i j = (b j : ℝ)})) := by ext p simp [S, IsFiniteTransportPlan] rw [heq] refine (isClosed_iInter fun i => isClosed_iInter fun j => ?_).inter ((isClosed_iInter fun i => ?_).inter (isClosed_iInter fun j => ?_)) · exact isClosed_le continuous_const (by fun_prop) · exact isClosed_eq (by fun_prop) continuous_const · exact isClosed_eq (by fun_prop) continuous_const have hbox : S ⊆ Set.Icc (0 : ι → κ → ℝ) (fun i j => (a i : ℝ)) := by intro p hp refine ⟨hp.1, ?_⟩ intro i j calc p i j ≤ ∑ k, p i k := Finset.single_le_sum (fun k _ => hp.1 i k) (Finset.mem_univ j) _ = (a i : ℝ) := hp.2.1 i have hcompact : IsCompact S := isCompact_Icc.of_isClosed_subset hclosed hbox obtain ⟨p, hp, hmin⟩ := hcompact.exists_isMinOn hne (show ContinuousOn (fun p : ι → κ → ℝ => ∑ i, ∑ j, c i j * p i j) S from (by fun_prop : Continuous (fun p : ι → κ → ℝ => ∑ i, ∑ j, c i j * p i j)).continuousOn) exact ⟨p, hp, fun q hq => hmin hq⟩ /-- Given [a cost matrix](hyp:c), [a feasible transport plan](hyp:hp) that [has the least total cost among all transport plans](hyp:hmin), and [any positive tolerance](hyp:hε), [there are row and column potentials whose sum never exceeds the cost of the corresponding entry and whose mass-weighted total comes within that tolerance of the plan's cost](goal). This is approximate linear-programming strong duality for the transportation problem, obtained from Causalean's conic duality theorem. -/ theorem finiteTransport_dual_certificate (a : ι → ℝ≥0) (b : κ → ℝ≥0) (c : ι → κ → ℝ) (p : ι → κ → ℝ) (hp : IsFiniteTransportPlan a b p) (hmin : ∀ q, IsFiniteTransportPlan a b q → (∑ i, ∑ j, c i j * p i j) ≤ ∑ i, ∑ j, c i j * q i j) {ε : ℝ} (hε : 0 < ε) : ∃ u : ι → ℝ, ∃ v : κ → ℝ, (∀ i j, u i + v j ≤ c i j) ∧ (∑ i, ∑ j, c i j * p i j) < (∑ i, (a i : ℝ) * u i) + (∑ j, (b j : ℝ) * v j) + ε := by -- Encode p in EuclideanSpace ℝ (ι × κ), with the nonnegative cone; -- constraints live in EuclideanSpace ℝ (ι ⊕ κ). Redundant total-mass -- equations preclude a naive Slater argument. Causalean's augmented-image -- qualification works: in a convergent augmented-image sequence, row sums -- bound all nonnegative entries, providing a convergent subsequence. -- Apply strong_duality_zero_gap and the defining property of a supremum. classical let K : ProperCone ℝ (EuclideanSpace ℝ (ι × κ)) := { carrier := {x | ∀ ij, 0 ≤ x ij} zero_mem' := by simp add_mem' := by intro x y hx hy ij exact add_nonneg (hx ij) (hy ij) smul_mem' := by intro r x hx ij exact mul_nonneg r.property (hx ij) isClosed' := by simp only [Set.ofPred_forall] exact isClosed_iInter fun ij => isClosed_le continuous_const (show Continuous (fun x : EuclideanSpace ℝ (ι × κ) => x ij) from by fun_prop) } let L : EuclideanSpace ℝ (ι × κ) →ₗ[ℝ] (EuclideanSpace ℝ (ι ⊕ κ)) := { toFun := fun x => WithLp.toLp 2 (Sum.elim (fun i => ∑ j, x (i, j)) (fun j => ∑ i, x (i, j))) map_add' := by intro x y apply PiLp.ext intro ij cases ij <;> simp [PiLp.add_apply, Finset.sum_add_distrib] map_smul' := by intro r x apply PiLp.ext intro ij cases ij <;> simp [PiLp.smul_apply, Finset.mul_sum] } let P : Causalean.PartialID.ConicProgram (EuclideanSpace ℝ (ι × κ)) (EuclideanSpace ℝ (ι ⊕ κ)) := { K := K A := L.toContinuousLinearMap b := WithLp.toLp 2 (Sum.elim (fun i => (a i : ℝ)) (fun j => (b j : ℝ))) c := WithLp.toLp 2 (fun ij => c ij.1 ij.2) } have hcost (x : EuclideanSpace ℝ (ι × κ)) : ⟪P.c, x⟫ = ∑ i, ∑ j, c i j * x (i, j) := by change inner ℝ (WithLp.toLp 2 (fun ij : ι × κ => c ij.1 ij.2)) x = _ rw [EuclideanSpace.inner_eq_star_dotProduct] simp [dotProduct, Fintype.sum_prod_type, mul_comm] have hfeas (x : EuclideanSpace ℝ (ι × κ)) : P.PrimalFeasible x ↔ IsFiniteTransportPlan a b (fun i j => x (i, j)) := by constructor · rintro ⟨heq, hK⟩ refine ⟨fun i j => hK (i, j), ?_, ?_⟩ · intro i exact congrArg (fun y : EuclideanSpace ℝ (ι ⊕ κ) => y (Sum.inl i)) heq · intro j exact congrArg (fun y : EuclideanSpace ℝ (ι ⊕ κ) => y (Sum.inr j)) heq · rintro ⟨hnonneg, hrow, hcol⟩ refine ⟨?_, fun ij => hnonneg ij.1 ij.2⟩ apply PiLp.ext intro ij cases ij with | inl i => exact hrow i | inr j => exact hcol j let x₀ : EuclideanSpace ℝ (ι × κ) := WithLp.toLp 2 (fun ij => p ij.1 ij.2) have hx₀ : P.PrimalFeasible x₀ := (hfeas x₀).2 hp have hlower (x : EuclideanSpace ℝ (ι × κ)) (hx : P.PrimalFeasible x) : ⟪P.c, x₀⟫ ≤ ⟪P.c, x⟫ := by simpa only [hcost] using hmin (fun i j => x (i, j)) ((hfeas x).1 hx) have hbdd : BddBelow ((fun x => ⟪P.c, x⟫) '' {x | P.PrimalFeasible x}) := by refine ⟨⟪P.c, x₀⟫, ?_⟩ rintro _ ⟨x, hx, rfl⟩ exact hlower x hx have hCQ : IsClosed P.augmentedImage := by apply IsSeqClosed.isClosed intro z zlim hz hlim choose x hx hxz using hz have htotal : Tendsto (fun n => ∑ i, (z n).1 (Sum.inl i)) atTop (𝓝 (∑ i, zlim.1 (Sum.inl i))) := by exact (show Continuous (fun z : EuclideanSpace ℝ (ι ⊕ κ) × ℝ => ∑ i, z.1 (Sum.inl i)) from by fun_prop).continuousAt.tendsto.comp hlim obtain ⟨M, hM⟩ := htotal.isCompact_insert_range.bddAbove have hbound (n : ℕ) (ij : ι × κ) : 0 ≤ x n ij ∧ x n ij ≤ M := by refine ⟨hx n ij, ?_⟩ calc x n ij ≤ ∑ j, x n (ij.1, j) := Finset.single_le_sum (fun j _ => hx n (ij.1, j)) (Finset.mem_univ ij.2) _ ≤ ∑ i, ∑ j, x n (i, j) := Finset.single_le_sum (fun i _ => Finset.sum_nonneg (fun j _ => hx n (i, j))) (Finset.mem_univ ij.1) _ = ∑ i, (z n).1 (Sum.inl i) := by have heq := congrArg Prod.fst (hxz n) simp only [P, L, LinearMap.toContinuousLinearMap] at heq rw [← heq] rfl _ ≤ M := hM (Set.mem_insert_of_mem _ (Set.mem_range_self n)) let B : Set (ι × κ → ℝ) := Set.Icc 0 (fun _ => M) have hxB : ∀ n, (x n).ofLp ∈ B := by intro n exact ⟨fun ij => (hbound n ij).1, fun ij => (hbound n ij).2⟩ obtain ⟨q, hq, φ, hφ, hqconv⟩ := (show IsCompact B from isCompact_Icc).tendsto_subseq hxB let xlim : EuclideanSpace ℝ (ι × κ) := WithLp.toLp 2 q have hxconv : Tendsto (x ∘ φ) atTop (𝓝 xlim) := (EuclideanSpace.equiv (ι × κ) ℝ).symm.continuous.continuousAt.tendsto.comp hqconv have hxlim : xlim ∈ P.K := by intro ij exact hq.1 ij have haug : Continuous (fun x : EuclideanSpace ℝ (ι × κ) => (P.A x, ⟪P.c, x⟫)) := by fun_prop have heqlim : (P.A xlim, ⟪P.c, xlim⟫) = zlim := by apply tendsto_nhds_unique (haug.continuousAt.tendsto.comp hxconv) have ht := hlim.comp hφ.tendsto_atTop simpa only [Function.comp_def, ← hxz] using ht exact ⟨xlim, hxlim, heqlim⟩ have hpair (x : EuclideanSpace ℝ (ι × κ)) (y : EuclideanSpace ℝ (ι ⊕ κ)) : ⟪P.A x, y⟫ = ∑ i, ∑ j, x (i, j) * (y (Sum.inl i) + y (Sum.inr j)) := by change inner ℝ (WithLp.toLp 2 (Sum.elim (fun i => ∑ j, x (i, j)) (fun j => ∑ i, x (i, j)))) y = _ rw [EuclideanSpace.inner_eq_star_dotProduct] simp only [dotProduct, Fintype.sum_sum_type, Pi.star_apply, star_trivial, Sum.elim_inl, Sum.elim_inr, Finset.mul_sum, mul_add, Finset.sum_add_distrib] rw [Finset.sum_comm (f := fun j i => y (Sum.inr j) * x (i, j))] simp only [mul_comm] have hdual (y : EuclideanSpace ℝ (ι ⊕ κ)) : P.DualFeasible y ↔ ∀ i j, y (Sum.inl i) + y (Sum.inr j) ≤ c i j := by constructor · intro hy i j let e : EuclideanSpace ℝ (ι × κ) := EuclideanSpace.single (i, j) 1 have he : e ∈ P.K := by intro ij simp only [e, PiLp.single_apply] split <;> norm_num have hred := ProperCone.mem_innerDual.mp hy he rw [inner_sub_right, ContinuousLinearMap.adjoint_inner_right] at hred have hecost : ⟪e, P.c⟫ = c i j := by change inner ℝ (EuclideanSpace.single (i, j) 1) (WithLp.toLp 2 _) = _ simp [EuclideanSpace.inner_single_left] have herow : ⟪P.A e, y⟫ = y (Sum.inl i) + y (Sum.inr j) := by rw [hpair] simp [e, PiLp.single_apply, Prod.ext_iff, ite_and] rw [hecost, herow] at hred linarith · intro hy apply ProperCone.mem_innerDual.mpr intro x hx rw [inner_sub_right, ContinuousLinearMap.adjoint_inner_right] rw [← real_inner_comm x P.c] rw [hcost, hpair] have heq : (∑ i, ∑ j, c i j * x (i, j)) - (∑ i, ∑ j, x (i, j) * (y (Sum.inl i) + y (Sum.inr j))) = ∑ i, ∑ j, (c i j - (y (Sum.inl i) + y (Sum.inr j))) * x (i, j) := by simp [sub_mul, Finset.sum_sub_distrib, mul_comm] rw [heq] exact Finset.sum_nonneg fun i _ => Finset.sum_nonneg fun j _ => mul_nonneg (sub_nonneg.mpr (hy i j)) (hx (i, j)) have hdual_ne : {y | P.DualFeasible y}.Nonempty := by let m := ∑ i, ∑ j, |c i j| refine ⟨WithLp.toLp 2 (Sum.elim (fun _ => -m) (fun _ => 0)), (hdual _).2 ?_⟩ intro i j change -m + 0 ≤ c i j have hcm : |c i j| ≤ m := by exact (Finset.single_le_sum (fun k _ => abs_nonneg (c i k)) (Finset.mem_univ j)).trans (Finset.single_le_sum (fun k _ => Finset.sum_nonneg (fun l _ => abs_nonneg (c k l))) (Finset.mem_univ i)) have hc := neg_abs_le (c i j) linarith have hvalue : P.primalValue = ⟪P.c, x₀⟫ := by apply le_antisymm · exact csInf_le hbdd ⟨x₀, hx₀, rfl⟩ · exact le_csInf ⟨⟪P.c, x₀⟫, x₀, hx₀, rfl⟩ (by rintro _ ⟨x, hx, rfl⟩ exact hlower x hx) have hgap := P.strong_duality_zero_gap ⟨x₀, hx₀⟩ hbdd hCQ have hlt : ⟪P.c, x₀⟫ - ε < P.dualValue := by rw [← hgap, hvalue]; linarith obtain ⟨t, ⟨y, hy, rfl⟩, hyt⟩ := exists_lt_of_lt_csSup (hdual_ne.image _) hlt refine ⟨fun i => y (Sum.inl i), fun j => y (Sum.inr j), ?_, ?_⟩ · exact (hdual y).1 hy · have hycost : ⟪P.b, y⟫ = (∑ i, (a i : ℝ) * y (Sum.inl i)) + ∑ j, (b j : ℝ) * y (Sum.inr j) := by change inner ℝ (WithLp.toLp 2 (Sum.elim (fun i => (a i : ℝ)) (fun j => (b j : ℝ)))) y = _ rw [EuclideanSpace.inner_eq_star_dotProduct] simp [dotProduct, Fintype.sum_sum_type, mul_comm] simp only [hycost, hcost] at hyt exact (sub_lt_iff_lt_add).1 hyt end Toolkit.Foundations.Optimization.Convex