getting tolerances right

src/functions/indHalfspace.jl

@@ -15,16 +15,13 @@ S = \\{x : \\langle a,x \\rangle \\leq b \\}.
 struct IndHalfspace{R <: Real, T <: AbstractArray{R}} <: ProximableFunction
   a::T
   b::R
+  norm_a::R
   function IndHalfspace{R, T}(a::T, b::R) where {R <: Real, T <: AbstractArray{R}}
     norm_a = norm(a)
     if norm_a == 0 && b < 0
         error("function is improper")
     end
-    if norm_a > 0
-        new(a/norm_a, b/norm_a)
-    else
-        new(a, b)
-    end
+    new(a, b, norm_a)
   end
 end
 
@@ -35,18 +32,18 @@ is_set(f::IndHalfspace) = true
 is_cone(f::IndHalfspace) = f.b == 0 || f.b == Inf
 
 function (f::IndHalfspace{R})(x::AbstractArray{R}) where R
-  if dot(f.a, x) <= f.b
+  if dot(f.a, x) - f.b <= eps(R)*f.norm_a*(1 + abs(f.b))
     return zero(R)
   end
   return R(Inf)
 end
 
 function prox!(y::AbstractArray{R}, f::IndHalfspace{R}, x::AbstractArray{R}, gamma::R=one(R)) where R
-  s = dot(f.a, x) - f.b
-  if s <= 0
-    copyto!(y, x)
+  s = dot(f.a, x)
+  if s > f.b
+    y .= x .- ((s - f.b)/f.norm_a^2) .* f.a
   else
-    y .= x .- s .* f.a
+    copyto!(y, x)
   end
   return zero(R)
 end
@@ -63,5 +60,5 @@ function prox_naive(f::IndHalfspace{R}, x::AbstractArray{R}, gamma::R=one(R)) wh
   if s <= 0
     return x, 0.0
   end
-  return x - s*f.a, 0.0
+  return x - (s/norm(f.a)^2)*f.a, 0.0
 end

src/functions/indHyperslab.jl

@@ -22,11 +22,7 @@ struct IndHyperslab{R <: Real, T <: AbstractArray{R}} <: ProximableFunction
     if (norm_a == 0 && (upp < 0 || low > 0)) || upp < low
         error("function is improper")
     end
-    if norm_a > 0
-        new(low/norm_a, a/norm_a, upp/norm_a, R(1.0))
-    else
-        new(low, a, upp, norm_a)
-    end
+    new(low, a, upp, norm_a)
   end
 end
 
@@ -43,19 +39,19 @@ is_cone(f::IndHyperslab{R}) where R =
 function (f::IndHyperslab{R})(x::AbstractArray{R}) where R
   s = dot(f.a, x)
   # if f.low <= s <= f.upp
-  tol = (R(1) + abs(s))*eps(R)
-  if f.low - s <= tol && s - f.upp <= tol
+  # tol = (R(1) + abs(s))*eps(R)
+  if f.low - s <= eps(R)*f.norm_a*(1 + abs(f.low)) && s - f.upp <= eps(R)*f.norm_a*(1 + abs(f.upp))
     return zero(R)
   end
   return R(Inf)
 end
 
 function prox!(y::AbstractArray{R}, f::IndHyperslab{R}, x::AbstractArray{R}, gamma::R=one(R)) where R
   s = dot(f.a, x)
-  if s < f.low
-    y .= x .- (s - f.low) .* f.a
-  elseif s > f.upp
-    y .= x .- (s - f.upp) .* f.a
+  if s < f.low && f.norm_a > 0
+    y .= x .- ((s - f.low)/f.norm_a^2) .* f.a
+  elseif s > f.upp && f.norm_a > 0
+    y .= x .- ((s - f.upp)/f.norm_a^2) .* f.a
   else
     copyto!(y, x)
   end
@@ -64,10 +60,10 @@ end
 
 function prox_naive(f::IndHyperslab{R}, x::AbstractArray{R}, gamma::R=one(R)) where R
   s = dot(f.a, x)
-  if s < f.low
-    return x - (s - f.low)*f.a, 0
-  elseif s > f.upp
-    return x - (s - f.upp)*f.a, 0
+  if s < f.low && f.norm_a > 0
+    return x - ((s - f.low)/norm(f.a)^2) * f.a, 0
+  elseif s > f.upp && f.norm_a > 0
+    return x - ((s - f.upp)/norm(f.a)^2) * f.a, 0
   else
     return x, 0
   end