diff --git a/Project.toml b/Project.toml
index bac4a3e..45b1a10 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,6 +1,6 @@
 name = "ProximalOperators"
 uuid = "a725b495-10eb-56fe-b38b-717eba820537"
-version = "0.15.3"
+version = "0.16.0"
 
 [deps]
 IterativeSolvers = "42fd0dbc-a981-5370-80f2-aaf504508153"
diff --git a/README.md b/README.md
index a6139b6..1f8d135 100644
--- a/README.md
+++ b/README.md
@@ -3,6 +3,7 @@
 [![Build status](https://github.com/JuliaFirstOrder/ProximalOperators.jl/workflows/CI/badge.svg?branch=master)](https://github.com/JuliaFirstOrder/ProximalOperators.jl/actions?query=workflow%3ACI+branch%3Amaster)
 [![codecov](https://codecov.io/gh/JuliaFirstOrder/ProximalOperators.jl/branch/master/graph/badge.svg)](https://codecov.io/gh/JuliaFirstOrder/ProximalOperators.jl)
 [![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.4020558.svg)](https://doi.org/10.5281/zenodo.4020558)
+[![Aqua QA](https://raw.githubusercontent.com/JuliaTesting/Aqua.jl/master/badge.svg)](https://github.com/JuliaTesting/Aqua.jl)
 
 [![](https://img.shields.io/badge/docs-stable-blue.svg)](https://juliafirstorder.github.io/ProximalOperators.jl/stable)
 [![](https://img.shields.io/badge/docs-latest-blue.svg)](https://juliafirstorder.github.io/ProximalOperators.jl/latest)
diff --git a/docs/Project.toml b/docs/Project.toml
index d8a086b..933768a 100644
--- a/docs/Project.toml
+++ b/docs/Project.toml
@@ -1,7 +1,7 @@
 [deps]
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
-ProximalOperators = "a725b495-10eb-56fe-b38b-717eba820537"
 ProximalCore = "dc4f5ac2-75d1-4f31-931e-60435d74994b"
+ProximalOperators = "a725b495-10eb-56fe-b38b-717eba820537"
 
 [compat]
 Documenter = "0.27"
diff --git a/src/ProximalOperators.jl b/src/ProximalOperators.jl
index 8142dbf..7b9575a 100644
--- a/src/ProximalOperators.jl
+++ b/src/ProximalOperators.jl
@@ -19,7 +19,6 @@ export prox, prox!, gradient, gradient!
 # Utilities
 
 include("utilities/approx_inequality.jl")
-include("utilities/tuples.jl")
 include("utilities/linops.jl")
 include("utilities/symmetricpacked.jl")
 include("utilities/uniformarrays.jl")
diff --git a/src/calculus/separableSum.jl b/src/calculus/separableSum.jl
index df3b48b..5142cf4 100644
--- a/src/calculus/separableSum.jl
+++ b/src/calculus/separableSum.jl
@@ -41,15 +41,23 @@ component_types(::Type{SeparableSum{T}}) where T = fieldtypes(T)
 
 (g::SeparableSum)(xs::Tuple) = sum(f(x) for (f, x) in zip(g.fs, xs))
 
-prox!(ys::Tuple, fs::Tuple, xs::Tuple, gamma::Number) = sum(prox!(y, f, x, gamma) for (y, f, x) in zip(ys, fs, xs))
+prox!(ys::Tuple, g::SeparableSum, xs::Tuple, gamma::Number) = sum(prox!(y, f, x, gamma) for (y, f, x) in zip(ys, g.fs, xs))
 
-prox!(ys::Tuple, fs::Tuple, xs::Tuple, gammas::Tuple) = sum(prox!(y, f, x, gamma) for (y, f, x, gamma) in zip(ys, fs, xs, gammas))
+prox!(ys::Tuple, g::SeparableSum, xs::Tuple, gammas::Tuple) = sum(prox!(y, f, x, gamma) for (y, f, x, gamma) in zip(ys, g.fs, xs, gammas))
 
-prox!(ys::Tuple, g::SeparableSum, xs::Tuple, gamma) = prox!(ys, g.fs, xs, gamma)
+function prox(g::SeparableSum, xs::Tuple, gamma=1)
+    ys = similar.(xs)
+    fys = prox!(ys, g, xs, gamma)
+    return ys, fys
+end
 
-gradient!(grads::Tuple, fs::Tuple, xs::Tuple) = sum(gradient!(grad, f, x) for (grad, f, x) in zip(grads, fs, xs))
+gradient!(grads::Tuple, g::SeparableSum, xs::Tuple) = sum(gradient!(grad, f, x) for (grad, f, x) in zip(grads, g.fs, xs))
 
-gradient!(grads::Tuple, g::SeparableSum, xs::Tuple) = gradient!(grads, g.fs, xs)
+function gradient(g::SeparableSum, xs::Tuple)
+    ys = similar.(xs)
+    fxs = gradient!(ys, g, xs)
+    return ys, fxs
+end
 
 function prox_naive(f::SeparableSum, xs::Tuple, gamma)
     fys = real(eltype(xs[1]))(0)
diff --git a/src/utilities/tuples.jl b/src/utilities/tuples.jl
deleted file mode 100644
index 166bcaa..0000000
--- a/src/utilities/tuples.jl
+++ /dev/null
@@ -1,9 +0,0 @@
-# Extension of some common Array operations to Tuple objects
-
-import Base: isapprox, similar, zero
-
-isapprox(x::Tuple, y::Tuple, args...; kwargs...) = all(isapprox.(x, y, args...; kwargs...))
-
-similar(x::Tuple) = similar.(x)
-
-zero(x::Tuple) = zero.(x)
diff --git a/test/Project.toml b/test/Project.toml
index a586ca3..1b3a6c1 100644
--- a/test/Project.toml
+++ b/test/Project.toml
@@ -1,5 +1,6 @@
 [deps]
-Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
-Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
+Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
+Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
diff --git a/test/runtests.jl b/test/runtests.jl
index 9043647..f91f581 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -18,6 +18,8 @@ using ProximalOperators:
     is_positively_homogeneous,
     is_support
 
+using Aqua
+
 function call_test(f, x::ArrayOrTuple{R}) where R <: Real
     try
         fx = @inferred f(x)
@@ -30,6 +32,8 @@ function call_test(f, x::ArrayOrTuple{R}) where R <: Real
     end
 end
 
+Base.zero(xs::Tuple) = Base.zero.(xs)
+
 # tests equality of the results of prox, prox! and prox_naive
 function prox_test(f, x::ArrayOrTuple{R}, gamma=1) where R <: Real
     y, fy = @inferred prox(f, x, gamma)
@@ -109,8 +113,8 @@ function predicates_test(f)
     @test !is_strongly_convex(f) || is_convex(f)
 end
 
-@testset "Sanity checks" begin
-    @test isempty(detect_unbound_args(ProximalOperators))
+@testset "Aqua" begin
+    Aqua.test_all(ProximalOperators; ambiguities=false)
 end
 
 @testset "Utilities" begin