diff --git a/src/ShiftedProximalOperators.jl b/src/ShiftedProximalOperators.jl
index 7db1fdf9..84a5813a 100644
--- a/src/ShiftedProximalOperators.jl
+++ b/src/ShiftedProximalOperators.jl
@@ -36,6 +36,7 @@ include("shiftedNormL0.jl")
 include("shiftedNormL0Box.jl")
 include("shiftedRootNormLhalf.jl")
 include("shiftedNormL1.jl")
+include("shiftedNormLinf.jl")
 include("shiftedGroupNormL2.jl")
 
 include("shiftedNormL1B2.jl")
@@ -45,6 +46,7 @@ include("shiftedIndBallL0BInf.jl")
 include("shiftedRootNormLhalfBox.jl")
 include("shiftedGroupNormL2Binf.jl")
 include("shiftedGroupNormL2Box.jl")
+include("shiftedNormLinfBox.jl")
 include("shiftedRank.jl")
 include("shiftedCappedl1.jl")
 include("shiftedNuclearnorm.jl")
diff --git a/src/shiftedNormLinf.jl b/src/shiftedNormLinf.jl
new file mode 100644
index 00000000..d18db1b1
--- /dev/null
+++ b/src/shiftedNormLinf.jl
@@ -0,0 +1,80 @@
+export ShiftedNormLinf
+
+"""
+Allows to compute the ‖.‖∞ operator with variable bounds: t ↦ λ ‖xk + sj + t‖∞
+
+The proximal operator is also provided. To do so, we use the algorithm proposed in Efficient Projections onto the ℓ1-Ball for Learning in High Dimensions¹
+
+¹https://ai.stanford.edu/~jduchi/projects/jd_ss_ys_l1.pdf
+"""
+
+mutable struct ShiftedNormLinf{
+  R <: Real,
+  V0 <: AbstractVector{R},
+  V1 <: AbstractVector{R},
+  V2 <: AbstractVector{R},
+} <: ShiftedProximableFunction
+  h::Conjugate{IndBallL1{R}}
+  xk::V0
+  sj::V1
+  sol::V2
+  shifted_twice::Bool
+  xsy::V2
+
+  function ShiftedNormLinf(
+    h::Conjugate{IndBallL1{R}},
+    xk::AbstractVector{R},
+    sj::AbstractVector{R},
+    shifted_twice::Bool,
+  ) where {R <: Real}
+    sol = similar(sj)
+    xsy = similar(sj)
+    new{R, typeof(xk), typeof(sj), typeof(sol)}(h, xk, sj, sol, shifted_twice, xsy)
+  end
+end
+
+shifted(h::Conjugate{IndBallL1{R}}, xk::AbstractVector{R}) where {R <: Real} = 
+    ShiftedNormLinf(h, xk, zero(xk), false)
+shifted(
+  ψ::ShiftedNormLinf{R, V0, V1, V2},
+  sj::AbstractVector{R},
+) where {R <: Real, V0 <: AbstractVector{R}, V1 <: AbstractVector{R}, V2 <: AbstractVector{R}} =
+  ShiftedNormLinf(ψ.h, ψ.xk, sj, true)
+
+fun_name(ψ::ShiftedNormLinf) = "shifted L∞ norm"
+fun_expr(ψ::ShiftedNormLinf) = "t ↦ λ ‖xk + sj + t‖∞"
+fun_params(ψ::ShiftedNormLinf) = "xk = $(ψ.xk)\n" * " "^14 * "sj = $(ψ.sj)\n" * " "^14
+
+
+function prox!(
+  y::AbstractVector{R},
+  ψ::ShiftedNormLinf{R, V0, V1, V2},
+  q::AbstractVector{R},
+  σ::R,
+) where {R <: Real, V0 <: AbstractVector{R}, V1 <: AbstractVector{R}, V2 <: AbstractVector{R}}
+  λ = ψ.h.f.r
+  @. ψ.sol = q + ψ.xk + ψ.sj
+
+  r = σ * λ
+  y .= ψ.sol .- _proj_l1ball(ψ.sol, r)
+  @. y -= (ψ.xk + ψ.sj)
+  val = zero(R)
+  @inbounds for i ∈ eachindex(y)
+    val = max(val, λ * abs(y[i] + ψ.xk[i] + ψ.sj[i]))
+  end
+  return val
+end
+
+function _proj_l1ball(v::AbstractVector{R}, r::R) where {R <: Real}
+  # Implements algorithm proposed in:
+  # Duchi et al. "Efficient Projections onto the ℓ₁-ball for Learning in High Dimensions",
+  if norm(v, 1) ≤ r
+    return copy(v)
+  end
+  μ = sort(abs.(v), rev = true)
+  cssμ = cumsum(μ)
+  list_inx = collect(1:length(μ))
+  rho = findlast(μ .* list_inx - (cssμ .- r) .+ eps(R) .> 0)
+  θ = (cssμ[rho] - r) / rho
+  return sign.(v) .* max.(abs.(v) .- θ, zero(R))
+end
\ No newline at end of file
diff --git a/src/shiftedNormLinfBox.jl b/src/shiftedNormLinfBox.jl
new file mode 100644
index 00000000..68b6eb97
--- /dev/null
+++ b/src/shiftedNormLinfBox.jl
@@ -0,0 +1,128 @@
+export ShiftedNormLinfBox
+
+"""
+Allows to compute the ‖.‖∞ operator with variable bounds: t ↦ λ ‖xk + sj + t‖∞ + χ({sj + t .∈ [l,u]})
+
+The proximal operator is also provided. To do so, we use the algorithm proposed in Efficient Projections onto the ℓ1-Ball for Learning in High Dimensions¹
+
+¹https://ai.stanford.edu/~jduchi/projects/jd_ss_ys_l1.pdf
+"""
+
+mutable struct ShiftedNormLinfBox{
+  R <: Real,
+  V0 <: AbstractVector{R},
+  V1 <: AbstractVector{R},
+  V2 <: AbstractVector{R},
+  V3,
+  V4,
+  VI <: AbstractArray{<:Integer},
+} <: ShiftedProximableFunction
+  h::Conjugate{IndBallL1{R}}
+  xk::V0
+  sj::V1
+  sol::V2
+  l::V3
+  u::V4
+  shifted_twice::Bool
+  selected::VI
+  xsy::V2
+
+  function ShiftedNormLinfBox(
+    h::Conjugate{IndBallL1{R}},
+    xk::AbstractVector{R},
+    sj::AbstractVector{R},
+    l,
+    u,
+    shifted_twice::Bool,
+    selected::AbstractArray{T},
+  ) where {R <: Real, T <: Integer}
+    sol = similar(xk)
+    xsy = similar(xk, length(selected))
+    if any(l .> u)
+      error("Error: at least one lower bound is greater than the upper bound.")
+    end
+    new{R, typeof(xk), typeof(sj), typeof(sol), typeof(l), typeof(u), typeof(selected)}(
+      h,
+      xk,
+      sj,
+      sol,
+      l,
+      u,
+      shifted_twice,
+      selected,
+      xsy,
+    )
+  end
+end
+
+shifted(
+  h::Conjugate{IndBallL1{R}},
+  xk::AbstractVector{R},
+  l,
+  u,
+  selected::AbstractArray{T} = 1:length(xk),
+) where {R <: Real, T <: Integer} =
+  ShiftedNormLinfBox(h, xk, zero(xk), l, u, false, selected)
+
+shifted(
+  ψ::ShiftedNormLinfBox{R, V0, V1, V2},
+  sj::AbstractVector{R},
+) where {R <: Real, V0 <: AbstractVector{R}, V1 <: AbstractVector{R}, V2 <: AbstractVector{R}} =
+  ShiftedNormLinfBox(ψ.h, ψ.xk, sj, ψ.l, ψ.u, true, ψ.selected)
+shifted(
+  ψ::ShiftedNormLinfBox{R, V0, V1, V2},
+  sj::AbstractVector{R},
+  l,
+  u,
+  selected::AbstractArray{T} = 1:length(sj),
+) where {R <: Real, T <: Integer, V0 <: AbstractVector{R}, V1 <: AbstractVector{R}, V2 <: AbstractVector{R}} =
+  ShiftedNormL1Box(ψ.h, ψ.xk, sj, l, u, true, selected)
+
+function (ψ::ShiftedNormLinfBox)(y)
+  tmp = ψ.xk .+ ψ.sj .+ y
+  val = ψ.h(tmp)
+  ϵ = √eps(eltype(y))
+  for i ∈ eachindex(y)
+    lower = isa(ψ.l, Real) ? ψ.l : ψ.l[i]
+    upper = isa(ψ.u, Real) ? ψ.u : ψ.u[i]
+    if !(lower - ϵ ≤ ψ.sj[i] + y[i] ≤ upper + ϵ)
+      return Inf
+    end
+  end
+  return val
+end
+
+fun_name(ψ::ShiftedNormLinfBox) = "shifted L∞ norm with box indicator"
+fun_expr(ψ::ShiftedNormLinfBox) = "t ↦ λ ‖xk + sj + t‖∞ + χ({sj + t .∈ [l,u]})"
+fun_params(ψ::ShiftedNormLinfBox) =
+  "xk = $(ψ.xk)\n" * " "^14 * "sj = $(ψ.sj)\n" * " "^14 * "l = $(ψ.l)\n" * " "^14 * "u = $(ψ.u)"
+
+function prox!(
+  y::AbstractVector{R},
+  ψ::ShiftedNormLinfBox{R, V0, V1, V2},
+  q::AbstractVector{R},
+  σ::R,
+) where {R <: Real, V0 <: AbstractVector{R}, V1 <: AbstractVector{R}, V2 <: AbstractVector{R}}
+  λ = ψ.h.f.r
+  @. ψ.sol = q + ψ.xk + ψ.sj
+
+  r = σ * λ
+  y .= ψ.sol .- _proj_l1ball(ψ.sol, r)
+
+  for i ∈ eachindex(y)
+    li = isa(ψ.l, Real) ? ψ.l : ψ.l[i]
+    ui = isa(ψ.u, Real) ? ψ.u : ψ.u[i]
+    si = ψ.sj[i]
+    qi = q[i]
+    if i ∈ ψ.selected
+      y[i] = min(max(y[i], li - si), ui - si)
+    else
+      y[i] = prox_zero(qi, li - si, ui - si)
+    end
+  end
+  val = zero(R)
+  @inbounds for i ∈ ψ.selected
+    val = max(val, λ * abs(y[i] + ψ.xk[i] + ψ.sj[i]))
+  end
+  return val
+end
\ No newline at end of file
diff --git a/test/runtests.jl b/test/runtests.jl
index 9855f2c4..d4ec7f58 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -1222,6 +1222,81 @@ for (op, shifted_op) ∈ zip((:Nuclearnorm,), (:ShiftedNuclearnorm,))
   end
 end
 
+for (op, shifted_op) ∈ zip((:NormLinf,), (:ShiftedNormLinf,))
+  @testset "$shifted_op" begin
+    ShiftedOp = eval(shifted_op)
+    Op = eval(op)
+    h = NormLinf(1.5)
+    n = 5
+    x = zeros(n)
+    ψ = shifted(h, x)
+    @test typeof(ψ) <: ShiftedOp
+    @test all(ψ.sj .== 0)
+    @test all(ψ.xk .== x)
+    # test value function
+    @test ψ(zeros(n)) == h(x)
+    # test shift update
+    y = rand(n)
+    shift!(ψ, y)
+    @test all(ψ.sj .== 0)
+    @test all(ψ.xk .== y)
+    # shift a shifted operator
+    s = ones(n) ./ 2
+    φ = shifted(ψ, s)
+    @test all(φ.sj .== s)
+    @test all(φ.xk .== y)
+    # test prox
+    q = randn(n)
+    ν = rand()
+    yp = similar(q)
+    val = prox!(yp, φ, q, ν)
+    @test val ≈ φ(yp)
+    @test val ≈ φ.h.f.r * norm(yp .+ φ.xk .+ φ.sj, Inf)
+  end
+end
+
+for (op, shifted_op) ∈ zip((:NormLinf,), (:ShiftedNormLinfBox,))
+  @testset "$shifted_op" begin
+    ShiftedOp = eval(shifted_op)
+    h = NormLinf(1.5)
+    n = 5
+    x = zeros(n)
+    l = -ones(n)
+    u = ones(n)
+    ψ = shifted(h, x, l, u)
+    @test typeof(ψ) <: ShiftedOp
+    @test all(ψ.sj .== 0)
+    @test all(ψ.xk .== x)
+    @test ψ.l == l
+    @test ψ.u == u
+    @test ψ.selected == 1:n
+    # test value function in the box
+    @test ψ(zeros(n)) == h(x)
+    # test value function out of the box -> Inf
+    @test ψ(2 .* ones(n)) == Inf
+    # test shift update
+    y = rand(n)
+    shift!(ψ, y)
+    @test all(ψ.sj .== 0)
+    @test all(ψ.xk .== y)
+    # shift a shifted operator
+    s = ones(n) ./ 2
+    φ = shifted(ψ, s)
+    @test all(φ.sj .== s)
+    @test all(φ.xk .== y)
+    @test φ.l == l
+    @test φ.u == u
+    # test prox
+    q = randn(n)
+    ν = rand()
+    yp = similar(q)
+    val = prox!(yp, φ, q, ν)
+    @test val ≈ φ(yp)
+    @test all(l .- s .- sqrt(eps()) .≤ yp .≤ u .- s .+ sqrt(eps()))
+    @test val ≈ φ.h.f.r * norm(yp[φ.selected] .+ φ.xk[φ.selected] .+ φ.sj[φ.selected], Inf)
+  end
+end
+
 for (op, shifted_op) ∈ zip((:GroupNormL2,), (:ShiftedGroupNormL2Box,))
   @testset "$shifted_op" begin
     ShiftedOp = eval(shifted_op)
diff --git a/test/test_allocs.jl b/test/test_allocs.jl
index aec1e857..328b8445 100644
--- a/test/test_allocs.jl
+++ b/test/test_allocs.jl
@@ -149,21 +149,17 @@ end
   end
 end
 
-for (op, shifted_op) ∈ zip((:GroupNormL2,), (:ShiftedGroupNormL2,), (:ShiftedGroupNormL2Box,))
-  λ = [1.0, 2.0]
-  idx = [1:500, 501:1000]
-  h = GroupNormL2(λ, idx)
-  n = 1000
-  xk = rand(n)
-  y = rand(n)
-  l = -3.0 * ones(n)
-  u = 4.0 * ones(n)
-
-  ψ = shifted(h, xk, l, u)
-  @test @wrappedallocs(ψ(y)) == 0
-  @test @wrappedallocs(prox!(y, ψ, y, 1.0)) == 0
-
-  ω = shifted(ψ, rand(n))
-  @test @wrappedallocs(ω(y)) == 0
-  @test @wrappedallocs(prox!(y, ω, y, 1.0)) == 0
-end
+# for op ∈ (:shiftedNormLinf,)
+#   h = NormLinf(1.0)
+#   n = 1000
+#   xk = rand(n)
+#   y = rand(n)
+
+#   ψ = shifted(h, xk)
+#   @test @wrappedallocs(ψ(y)) == 0
+#   @test @wrappedallocs(prox!(y, ψ, y, 1.0)) == 0
+
+#   ψ = shifted(h, xk, -3.0 * ones(n), 4.0 * ones(n), rand(1:n, Int(n / 2)))
+#   @test @wrappedallocs(ψ(y)) == 0
+#   @test @wrappedallocs(prox!(y, ψ, y, 1.0)) == 0
+# end
diff --git a/test/testsbox.jl b/test/testsbox.jl
index 03cb1d12..2e066362 100644
--- a/test/testsbox.jl
+++ b/test/testsbox.jl
@@ -365,3 +365,76 @@ for (op, shifted_op) ∈ zip((:GroupNormL2,), (:ShiftedGroupNormL2Box,))
     @test ω.selected == selected
   end
 end
+
+for (op, shifted_op) ∈ zip((:NormLinf,), (:ShiftedNormLinfBox,))
+
+  function _proj_l1ball(v::AbstractVector{R}, r::R) where {R <: Real}
+    # Implements algorithm proposed in:
+    # Duchi et al. "Efficient Projections onto the ℓ₁-ball for Learning in High Dimensions",
+    if norm(v, 1) ≤ r
+      return copy(v)
+    end
+    μ = sort(abs.(v), rev = true)
+    cssμ = cumsum(μ)
+    list_inx = collect(1:length(μ))
+    rho = findlast(μ .* list_inx - (cssμ .- r) .+ eps(R) .> 0)
+    θ = (cssμ[rho] - r) / rho
+    return sign.(v) .* max.(abs.(v) .- θ, zero(R))
+  end
+
+  @testset "$shifted_op" begin
+    ## Testing ShiftedNormLinfBox
+    # Looking for argmin_t obj(t) = 1/(2σ) * ‖t-q‖² + λ * ‖x+s+t‖∞ + χ{s+t ∈ [l,u]}
+    σ = 1.0
+    λ = 1.5
+    l = [-1.0, -1.0]
+    u = [1.0, 1.0]
+    s = [0.0, 0.0]
+    x = [[0.0, 0.0], [0.0, 0.0], [0.0, 0.0], [0.5, 0.5]]
+    q = [[3.0, 4.0], [0.1, 0.1], [-3.0, -4.0], [2.5, 3.5]]
+    # sol = [[1.0, 1.0], [0.1, 0.1] .- _proj_l1ball([0.1, 0.1], σ * λ), [-1.0, -1.0], [1.0, 1.0]]
+    sol = [[1.0, 1.0], [0.0, 0.0], [-1.0, -1.0], [1.0, 1.0]]
+    # Case 1 : prox out of box -> clip u-s
+    # Case 2 : q small, prox in box -> no clip
+    # Case 3 : prox out of box -> clip l-s
+    # Case 4 : xk ≠ 0, prox out of box -> clip u-s
+    for i = 1:4
+      h = NormLinf(λ)
+      ψ = shifted(h, x[i], l, u)
+      ω = shifted(ψ, s)
+      y = similar(q[i])
+      val = prox!(y, ω, q[i], σ)
+      @test isapprox(y, sol[i], atol = 1.0e-8)
+      @test all(l .- s .- sqrt(eps()) .≤ y .≤ u .- s .+ sqrt(eps()))
+      @test isapprox(val, ω(y), atol = 1.0e-8)
+    end
+  end
+  @testset "ShiftedNormLinfBox with selected" begin
+    σ = 1.0
+    λ = 1.5
+    l = [-1.0, -1.0, -1.0]
+    u = [1.0, 1.0, 1.0]
+    s = [0.0, 0.0, 0.0]
+    x = [0.0, 0.0, 0.0]
+    selected = [1, 2]  # index 3 excluded from the norm
+    q   = [[3.0, 4.0, 0.0], [3.0, 4.0, 5.0], [3.0, 4.0, -5.0]]
+    sol = [[1.0, 1.0, 0.0], [1.0, 1.0, 1.0], [1.0, 1.0, -1.0]]
+    # Case 1 : q[3] = 0 ∈ [l-s,u-s] -> prox_zero(0,-1,1) = 0
+    # Case 2 : q[3] = 5 > u-s = 1   -> prox_zero(5,-1,1) = 1
+    # Case 3 : q[3] = -5 < l-s = -1 -> prox_zero(-5,-1,1) = -1
+    for i = 1:3
+      h = NormLinf(λ)
+      ψ = shifted(h, x, l, u, selected)
+      ω = shifted(ψ, s)
+      y = similar(q[i])
+      val = prox!(y, ω, q[i], σ)
+      @test isapprox(y, sol[i], atol = 1.0e-8)
+      @test isapprox(val, ω(y), atol = 1.0e-8)
+    end
+    # check selected is propagated after shift
+    h = NormLinf(λ)
+    ψ = shifted(h, x, l, u, selected)
+    ω = shifted(ψ, s)
+    @test ω.selected == selected
+  end
+end
\ No newline at end of file