diff --git a/Project.toml b/Project.toml
index 9c1a305c..898241a9 100644
--- a/Project.toml
+++ b/Project.toml
@@ -8,15 +8,16 @@ projects = ["docs", "test", "benchmark"]
 [deps]
 IterativeSolvers = "42fd0dbc-a981-5370-80f2-aaf504508153"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
-OSQP = "ab2f91bb-94b4-55e3-9ba0-7f65df51de79"
 ProximalCore = "dc4f5ac2-75d1-4f31-931e-60435d74994b"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 TSVD = "9449cd9e-2762-5aa3-a617-5413e99d722e"
 
 [weakdeps]
+OSQP = "ab2f91bb-94b4-55e3-9ba0-7f65df51de79"
 RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd"
 
 [extensions]
+ProximalOperatorsOSQPExt = "OSQP"
 RecursiveArrayToolsExt = "RecursiveArrayTools"
 
 [compat]
diff --git a/src/functions/indPolyhedralOSQP.jl b/ext/ProximalOperatorsOSQPExt.jl
similarity index 64%
rename from src/functions/indPolyhedralOSQP.jl
rename to ext/ProximalOperatorsOSQPExt.jl
index 6fb67d7f..3280d9eb 100644
--- a/src/functions/indPolyhedralOSQP.jl
+++ b/ext/ProximalOperatorsOSQPExt.jl
@@ -1,52 +1,43 @@
-# IndPolyhedral: OSQP implementation
+module ProximalOperatorsOSQPExt
 
+using LinearAlgebra
+using SparseArrays
 using OSQP
 
-struct IndPolyhedralOSQP{R} <: IndPolyhedral
-    l::AbstractVector{R}
-    A::AbstractMatrix{R}
-    u::AbstractVector{R}
-    mod::OSQP.Model
-    function IndPolyhedralOSQP{R}(
-        l::AbstractVector{R}, A::AbstractMatrix{R}, u::AbstractVector{R}
-    ) where R
-        m, n = size(A)
-        mod = OSQP.Model()
-        if !all(l .<= u)
-            error("function is improper (are some bounds inverted?)")
-        end
-        OSQP.setup!(mod; P=SparseMatrixCSC{R}(I, n, n), l=l, A=sparse(A), u=u, verbose=false,
-            eps_abs=eps(R), eps_rel=eps(R),
-            eps_prim_inf=eps(R), eps_dual_inf=eps(R))
-        new(l, A, u, mod)
-    end
-end
-
-# properties
-
-is_proximable(::Type{<:IndPolyhedralOSQP}) = false
+using ProximalOperators
+using ProximalOperators: IndPolyhedralOSQP
+import ProximalCore: prox, prox!
 
 # constructors
 
-IndPolyhedralOSQP(
+function ProximalOperators.IndPolyhedralOSQP(
     l::AbstractVector{R}, A::AbstractMatrix{R}, u::AbstractVector{R}
-) where R =
-    IndPolyhedralOSQP{R}(l, A, u)
+) where R
+    m, n = size(A)
+    if !all(l .<= u)
+        error("function is improper (are some bounds inverted?)")
+    end
+    mod = OSQP.Model()
+    OSQP.setup!(mod; P=SparseMatrixCSC{R}(I, n, n), l=l, A=sparse(A), u=u, verbose=false,
+        eps_abs=eps(R), eps_rel=eps(R),
+        eps_prim_inf=eps(R), eps_dual_inf=eps(R))
+    return IndPolyhedralOSQP{R, typeof(mod)}(l, A, u, mod)
+end
 
-IndPolyhedralOSQP(
+ProximalOperators.IndPolyhedralOSQP(
     l::AbstractVector{R}, A::AbstractMatrix{R}, u::AbstractVector{R},
     xmin::AbstractVector{R}, xmax::AbstractVector{R}
 ) where R =
     IndPolyhedralOSQP([l; xmin], [A; I], [u; xmax])
 
-IndPolyhedralOSQP(
+ProximalOperators.IndPolyhedralOSQP(
     l::AbstractVector{R}, A::AbstractMatrix{R}, args...
 ) where R =
     IndPolyhedralOSQP(
         l, SparseMatrixCSC(A), R(Inf).*ones(R, size(A, 1)), args...
     )
 
-IndPolyhedralOSQP(
+ProximalOperators.IndPolyhedralOSQP(
     A::AbstractMatrix{R}, u::AbstractVector{R}, args...
 ) where R =
     IndPolyhedralOSQP(
@@ -81,7 +72,7 @@ end
 # dual problem is: minimize_y (1/2)||-A'y||^2 - x'A'y + g*(y)
 # can solve with (fast) dual proximal gradient method
 
-function prox_naive(f::IndPolyhedralOSQP, x, gamma)
+function ProximalOperators.prox_naive(f::IndPolyhedralOSQP, x, gamma)
     R = eltype(x)
     y = zeros(R, size(f.A, 1)) # dual vector
     y1 = y
@@ -99,3 +90,5 @@ function prox_naive(f::IndPolyhedralOSQP, x, gamma)
     p = -f.A'*y + x
     return p, R(0)
 end
+
+end # module ProximalOperatorsOSQPExt
diff --git a/src/functions/indPolyhedral.jl b/src/functions/indPolyhedral.jl
index 8405c2c3..989db847 100644
--- a/src/functions/indPolyhedral.jl
+++ b/src/functions/indPolyhedral.jl
@@ -14,6 +14,9 @@ S = \\{ x : x_\\min \\leq x \\leq x_\\max, l \\leq Ax \\leq u \\}.
 ```
 Matrix `A` is a mandatory argument; when any of the bounds is not provided,
 it is assumed to be (plus or minus) infinity.
+
+The default `solver=:osqp` backend is provided by a package extension: load it
+with `using OSQP` before constructing the object, otherwise an error is raised.
 """
 function IndPolyhedral(args...; solver=:osqp)
     if solver == :osqp
@@ -23,6 +26,34 @@ function IndPolyhedral(args...; solver=:osqp)
     end
 end
 
-# including concrete types
+# IndPolyhedral: OSQP implementation
+#
+# The struct below is defined here so that the `solver=:osqp` dispatch above can
+# refer to it, but everything that actually needs OSQP -- the constructors,
+# `prox!`, the function evaluation and `prox_naive` -- lives in the package
+# extension `ext/ProximalOperatorsOSQPExt.jl`, loaded once OSQP is available
+# (`using OSQP`). Without OSQP loaded, constructing an `IndPolyhedralOSQP`
+# (directly or via `IndPolyhedral(...; solver=:osqp)`) raises an informative
+# error.
+
+struct IndPolyhedralOSQP{R, M} <: IndPolyhedral
+    l::AbstractVector{R}
+    A::AbstractMatrix{R}
+    u::AbstractVector{R}
+    mod::M
+    # Explicit inner constructor: suppresses the auto-generated 4-positional-arg
+    # outer constructor, which would otherwise shadow the `(l, A, xmin, xmax)`
+    # constructor added by ProximalOperatorsOSQPExt.
+    IndPolyhedralOSQP{R, M}(l, A, u, mod) where {R, M} = new{R, M}(l, A, u, mod)
+end
 
-include("indPolyhedralOSQP.jl")
+is_proximable(::Type{<:IndPolyhedralOSQP}) = false
+
+# The real constructors are added to this function by ProximalOperatorsOSQPExt;
+# this fallback only fires when OSQP is not loaded.
+function IndPolyhedralOSQP(args...; kwargs...)
+    error(
+        "IndPolyhedralOSQP requires the OSQP package: run `using OSQP` before " *
+        "constructing IndPolyhedralOSQP(...) or IndPolyhedral(...; solver=:osqp)."
+    )
+end
diff --git a/test/Project.toml b/test/Project.toml
index 07c19fe7..639790d4 100644
--- a/test/Project.toml
+++ b/test/Project.toml
@@ -2,6 +2,7 @@
 Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+OSQP = "ab2f91bb-94b4-55e3-9ba0-7f65df51de79"
 ProximalCore = "dc4f5ac2-75d1-4f31-931e-60435d74994b"
 ProximalOperators = "a725b495-10eb-56fe-b38b-717eba820537"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
diff --git a/test/test_indPolyhedral.jl b/test/test_indPolyhedral.jl
index 0d0a0bb0..4322f52d 100644
--- a/test/test_indPolyhedral.jl
+++ b/test/test_indPolyhedral.jl
@@ -1,4 +1,5 @@
 using ProximalOperators
+using OSQP  # loads ProximalOperatorsOSQPExt, enabling IndPolyhedral(...; solver=:osqp)
 using Test
 
 @testset "IndPolyhedral" begin