pointsource_algs: comparison src/subproblem/unconstrained.rs

-:aec67cdd6b14
+:efa60bc4f743
 use crate::types::*;
 use super::{
 InnerMethod,
 InnerSettings
 };
+use super::l1squared_nonneg::max_interval_dist_to_zero;
 /// Compute the proximal operator of $x \mapsto |x|$, i.e., the soft-thresholding operator.
 #[inline]
 #[replace_float_literals(F::cast_from(literal))]
-fn soft_thresholding<F : Float>(v : F, λ : F) -> F {
+pub(crate) fn soft_thresholding<F : Float>(v : F, λ : F) -> F {
 if v > λ {
 v - λ
 } else if v < -λ {
 v + λ
 } else {
 }
 /// Returns the ∞-norm minimal subdifferential of $x ↦  x^⊤Ax - g^⊤ x + λ\|x\|₁$ at $x$.
 ///
 /// `v` will be modified and cannot be trusted to contain useful values afterwards.
-#[replace_float_literals(F::cast_from(literal).to_nalgebra_mixed())]
+#[replace_float_literals(F::cast_from(literal))]
-fn min_subdifferential<F : Float + ToNalgebraRealField>(
+fn min_subdifferential<F : Float + nalgebra::RealField>(
-v : &mut DVector<F::MixedType>,
+v : &mut DVector<F>,
-mA : &DMatrix<F::MixedType>,
+mA : &DMatrix<F>,
-x : &DVector<F::MixedType>,
+x : &DVector<F>,
-g : &DVector<F::MixedType>,
+g : &DVector<F>,
-λ : F::MixedType
+λ : F
 ) -> F {
 v.copy_from(g);
 mA.gemv(v, 1.0, x, -1.0);   // v =  Ax - g
 let mut val = 0.0;
 for (&v_i, &x_i) in izip!(v.iter(), x.iter()) {
-// The subdifferential at x is $Ax - g + λ ∂‖·‖₁(x)$.
+let (mut lb, mut ub) = (v_i, v_i);
-val = val.max(match x_i.partial_cmp(&0.0) {
+match x_i.partial_cmp(&0.0) {
-Some(Greater) => v_i + λ,
+Some(Greater) => { lb += λ; ub += λ },
-Some(Less) => v_i - λ,
+Some(Less) => { lb -= λ; ub -= λ },
-Some(Equal) => soft_thresholding(v_i, λ),
+Some(Equal) => { lb -= λ; ub += λ },
-None => F::MixedType::nan(),
+None => {},
-})
+}
+val = max_interval_dist_to_zero(val, lb, ub);
 }
-F::from_nalgebra_mixed(val)
+val
 }
 /// Forward-backward splitting implementation of [`quadratic_unconstrained`].
 /// For detailed documentation of the inputs and outputs, refer to there.
 });
 iters +=1;
 backup.map(|_| {
-min_subdifferential(&mut v, mA, x, g, λ)
+F::from_nalgebra_mixed(min_subdifferential(&mut v, mA, x, g, λ))
 })
 });
 iters
 }
 iters += 1;
 // 4. Report solution quality
 state.if_verbose(|| {
 // Calculate subdifferential at the FB step `x` that hasn't yet had `s` yet added.
-min_subdifferential(&mut v, mA, x, g, λ)
+F::from_nalgebra_mixed(min_subdifferential(&mut v, mA, x, g, λ))
 })
 });
 res.map(|_| iters)
 }
 /// This function applies an iterative method for the solution of the problem
 /// <div>$$
-///     \min_{x ∈ ℝ^n} \frac{1}{2} x^⊤Ax - g^⊤ x + λ\|x\|₁ + c.
+///     \min_{x ∈ ℝ^n} \frac{1}{2} x^⊤Ax - g^⊤ x + λ\|x\|₁.
 /// $$</div>
 /// Semismooth Newton or forward-backward are supported based on the setting in `method`.
 /// The parameter `mA` is matrix $A$, and `g` and `λ` are as in the mathematical formulation.
 /// The constant $c$ does not need to be provided. The step length parameter is `τ` while
 /// `x` contains the initial iterate and on return the final one. The `iterator` controls
 /// stopping. The “verbose” value output by all methods is the $ℓ\_∞$ distance of some
 /// subdifferential of the objective to zero.
 ///
 /// This function returns the number of iterations taken.
 pub fn quadratic_unconstrained<F, I>(
-method : InnerMethod,
 mA : &DMatrix<F::MixedType>,
 g : &DVector<F::MixedType>,
-//c_ : F,
 λ : F,
 x : &mut DVector<F::MixedType>,
-τ : F,
+mA_normest : F,
+inner : &InnerSettings<F>,
 iterator : I
 ) -> usize
 where F : Float + ToNalgebraRealField,
 I : AlgIteratorFactory<F>
 {
+let inner_τ = inner.fb_τ0 / mA_normest;
-match method {
+match inner.method {
-InnerMethod::FB =>
+InnerMethod::FB | InnerMethod::Default =>
-quadratic_unconstrained_fb(mA, g, λ, x, τ, iterator),
+quadratic_unconstrained_fb(mA, g, λ, x, inner_τ, iterator),
 InnerMethod::SSN =>
-quadratic_unconstrained_ssn(mA, g, λ, x, τ, iterator).unwrap_or_else(|e| {
+quadratic_unconstrained_ssn(mA, g, λ, x, inner_τ, iterator).unwrap_or_else(|e| {
 println!("{}", format!("{e}. Using FB fallback.").red());
 let ins = InnerSettings::<F>::default();
-quadratic_unconstrained_fb(mA, g, λ, x, τ, ins.iterator_options)
+quadratic_unconstrained_fb(mA, g, λ, x, inner_τ, ins.iterator_options)
-})
+}),
+InnerMethod::PDPS => unimplemented!(),
 }
 }

Mercurial > repos > pointsource_algs / file comparison

comparison: src/subproblem/unconstrained.rs

src/subproblem/unconstrained.rs