diff --git a/docs/src/submodules/Nonlinear/reference.md b/docs/src/submodules/Nonlinear/reference.md
index 842db49713..468f44cdfe 100644
--- a/docs/src/submodules/Nonlinear/reference.md
+++ b/docs/src/submodules/Nonlinear/reference.md
@@ -72,6 +72,11 @@ Nonlinear.AbstractAutomaticDifferentiation
 Nonlinear.ExprGraphOnly
 Nonlinear.SparseReverseMode
 Nonlinear.SymbolicMode
+Nonlinear.QPBlockData
+Nonlinear.add_constraint_jacobian_product
+Nonlinear.add_constraint_jacobian_transpose_product
+Nonlinear.add_hessian_lagrangian_product
+
 ```
 
 ## Data-structure
diff --git a/src/Nonlinear/Nonlinear.jl b/src/Nonlinear/Nonlinear.jl
index 9a5e4bd8d9..e057052029 100644
--- a/src/Nonlinear/Nonlinear.jl
+++ b/src/Nonlinear/Nonlinear.jl
@@ -42,4 +42,6 @@ include("evaluator.jl")
 include("ReverseAD/ReverseAD.jl")
 include("SymbolicAD/SymbolicAD.jl")
 
+include("qp_block_data.jl")
+
 end  # module
diff --git a/src/Nonlinear/qp_block_data.jl b/src/Nonlinear/qp_block_data.jl
new file mode 100644
index 0000000000..85e6ad8ce4
--- /dev/null
+++ b/src/Nonlinear/qp_block_data.jl
@@ -0,0 +1,794 @@
+# Copyright (c) 2013: Iain Dunning, Miles Lubin, and contributors
+#
+# Use of this source code is governed by an MIT-style license that can be found
+# in the LICENSE.md file or at https://opensource.org/licenses/MIT.
+
+# This file is adapted from `Ipopt.jl/ext/IpoptMathOptInterfaceExt/utils.jl`.
+#
+# Unlike the Ipopt version, a variable is treated as a parameter if and only
+# if its index is a key of the `parameters` dictionary, instead of an
+# index-offset convention. Parameters must therefore be registered in
+# `parameters` before any structure query, but their values may be updated
+# freely between function evaluations.
+
+@enum(
+    _FunctionType,
+    _kFunctionTypeVariableIndex,
+    _kFunctionTypeScalarAffine,
+    _kFunctionTypeScalarQuadratic,
+)
+
+function _function_type_to_func(::Type{T}, k::_FunctionType) where {T}
+    if k == _kFunctionTypeVariableIndex
+        return MOI.VariableIndex
+    elseif k == _kFunctionTypeScalarAffine
+        return MOI.ScalarAffineFunction{T}
+    else
+        @assert k == _kFunctionTypeScalarQuadratic
+        return MOI.ScalarQuadraticFunction{T}
+    end
+end
+
+_function_info(::MOI.VariableIndex) = _kFunctionTypeVariableIndex
+_function_info(::MOI.ScalarAffineFunction) = _kFunctionTypeScalarAffine
+_function_info(::MOI.ScalarQuadraticFunction) = _kFunctionTypeScalarQuadratic
+
+@enum(
+    _BoundType,
+    _kBoundTypeLessThan,
+    _kBoundTypeGreaterThan,
+    _kBoundTypeEqualTo,
+    _kBoundTypeInterval,
+)
+
+_set_info(s::MOI.LessThan) = _kBoundTypeLessThan, -Inf, s.upper
+_set_info(s::MOI.GreaterThan) = _kBoundTypeGreaterThan, s.lower, Inf
+_set_info(s::MOI.EqualTo) = _kBoundTypeEqualTo, s.value, s.value
+_set_info(s::MOI.Interval) = _kBoundTypeInterval, s.lower, s.upper
+
+function _bound_type_to_set(::Type{T}, k::_BoundType) where {T}
+    if k == _kBoundTypeEqualTo
+        return MOI.EqualTo{T}
+    elseif k == _kBoundTypeLessThan
+        return MOI.LessThan{T}
+    elseif k == _kBoundTypeGreaterThan
+        return MOI.GreaterThan{T}
+    else
+        @assert k == _kBoundTypeInterval
+        return MOI.Interval{T}
+    end
+end
+
+"""
+    QPBlockData{T}()
+
+A data structure holding an affine or quadratic objective and a block of
+affine and quadratic constraints, together with methods to evaluate them
+following the [`MOI.AbstractNLPEvaluator`](@ref) callback conventions.
+
+This is a helper for solvers that pass affine and quadratic constraints to
+the solver through the same callbacks as an [`MOI.AbstractNLPEvaluator`](@ref)
+(for example, Ipopt and MadNLP).
+
+## Parameters
+
+A variable is treated as a parameter if and only if its index is a key of the
+`parameters` dictionary, which maps the raw `MOI.VariableIndex` value of the
+parameter to its current value. Register every parameter in `parameters`
+before querying any structure; the values may be updated freely between
+function evaluations.
+"""
+mutable struct QPBlockData{T}
+    objective::Union{MOI.ScalarAffineFunction{T},MOI.ScalarQuadraticFunction{T}}
+    objective_function_type::_FunctionType
+    constraints::Vector{
+        Union{MOI.ScalarAffineFunction{T},MOI.ScalarQuadraticFunction{T}},
+    }
+    g_L::Vector{T}
+    g_U::Vector{T}
+    mult_g::Vector{Union{Nothing,T}}
+    function_type::Vector{_FunctionType}
+    bound_type::Vector{_BoundType}
+    parameters::Dict{Int64,T}
+
+    function QPBlockData{T}() where {T}
+        return new(
+            zero(MOI.ScalarQuadraticFunction{T}),
+            _kFunctionTypeScalarAffine,
+            Union{MOI.ScalarAffineFunction{T},MOI.ScalarQuadraticFunction{T}}[],
+            T[],
+            T[],
+            Union{Nothing,T}[],
+            _FunctionType[],
+            _BoundType[],
+            Dict{Int64,T}(),
+        )
+    end
+end
+
+_is_parameter(v::MOI.VariableIndex, p::Dict) = haskey(p, v.value)
+
+function _value(v::MOI.VariableIndex, x, p::Dict)
+    return _is_parameter(v, p) ? p[v.value] : x[v.value]
+end
+
+function _eval_function(
+    f::MOI.ScalarQuadraticFunction{T},
+    x::AbstractVector{T},
+    p::Dict{Int64,T},
+)::T where {T}
+    y = f.constant
+    for term in f.affine_terms
+        y += term.coefficient * _value(term.variable, x, p)
+    end
+    for term in f.quadratic_terms
+        v1 = _value(term.variable_1, x, p)
+        v2 = _value(term.variable_2, x, p)
+        if term.variable_1 == term.variable_2
+            y += term.coefficient * v1 * v2 / 2
+        else
+            y += term.coefficient * v1 * v2
+        end
+    end
+    return y
+end
+
+function _eval_function(
+    f::MOI.ScalarAffineFunction{T},
+    x::AbstractVector{T},
+    p::Dict{Int64,T},
+)::T where {T}
+    y = f.constant
+    for term in f.terms
+        y += term.coefficient * _value(term.variable, x, p)
+    end
+    return y
+end
+
+function _eval_dense_gradient(
+    ∇f::AbstractVector{T},
+    f::MOI.ScalarQuadraticFunction{T},
+    x::AbstractVector{T},
+    p::Dict{Int64,T},
+)::Nothing where {T}
+    for term in f.affine_terms
+        if !_is_parameter(term.variable, p)
+            ∇f[term.variable.value] += term.coefficient
+        end
+    end
+    for term in f.quadratic_terms
+        if !_is_parameter(term.variable_1, p)
+            v = _value(term.variable_2, x, p)
+            ∇f[term.variable_1.value] += term.coefficient * v
+        end
+        if term.variable_1 != term.variable_2 &&
+           !_is_parameter(term.variable_2, p)
+            v = _value(term.variable_1, x, p)
+            ∇f[term.variable_2.value] += term.coefficient * v
+        end
+    end
+    return
+end
+
+function _eval_dense_gradient(
+    ∇f::AbstractVector{T},
+    f::MOI.ScalarAffineFunction{T},
+    x::AbstractVector{T},
+    p::Dict{Int64,T},
+)::Nothing where {T}
+    for term in f.terms
+        if !_is_parameter(term.variable, p)
+            ∇f[term.variable.value] += term.coefficient
+        end
+    end
+    return
+end
+
+function _append_sparse_gradient_structure!(
+    f::MOI.ScalarQuadraticFunction,
+    J,
+    row,
+    p::Dict,
+)
+    for term in f.affine_terms
+        if !_is_parameter(term.variable, p)
+            push!(J, (row, term.variable.value))
+        end
+    end
+    for term in f.quadratic_terms
+        if !_is_parameter(term.variable_1, p)
+            push!(J, (row, term.variable_1.value))
+        end
+        if term.variable_1 != term.variable_2 &&
+           !_is_parameter(term.variable_2, p)
+            push!(J, (row, term.variable_2.value))
+        end
+    end
+    return
+end
+
+function _append_sparse_gradient_structure!(
+    f::MOI.ScalarAffineFunction,
+    J,
+    row,
+    p::Dict,
+)
+    for term in f.terms
+        if !_is_parameter(term.variable, p)
+            push!(J, (row, term.variable.value))
+        end
+    end
+    return
+end
+
+function _eval_sparse_gradient(
+    ∇f::AbstractVector{T},
+    f::MOI.ScalarQuadraticFunction{T},
+    x::AbstractVector{T},
+    p::Dict{Int64,T},
+)::Int where {T}
+    i = 0
+    for term in f.affine_terms
+        if !_is_parameter(term.variable, p)
+            i += 1
+            ∇f[i] = term.coefficient
+        end
+    end
+    for term in f.quadratic_terms
+        if !_is_parameter(term.variable_1, p)
+            v = _value(term.variable_2, x, p)
+            i += 1
+            ∇f[i] = term.coefficient * v
+        end
+        if term.variable_1 != term.variable_2 &&
+           !_is_parameter(term.variable_2, p)
+            v = _value(term.variable_1, x, p)
+            i += 1
+            ∇f[i] = term.coefficient * v
+        end
+    end
+    return i
+end
+
+function _eval_sparse_gradient(
+    ∇f::AbstractVector{T},
+    f::MOI.ScalarAffineFunction{T},
+    x::AbstractVector{T},
+    p::Dict{Int64,T},
+)::Int where {T}
+    i = 0
+    for term in f.terms
+        if !_is_parameter(term.variable, p)
+            i += 1
+            ∇f[i] = term.coefficient
+        end
+    end
+    return i
+end
+
+function _append_sparse_hessian_structure!(
+    f::MOI.ScalarQuadraticFunction,
+    H,
+    p::Dict,
+)
+    for term in f.quadratic_terms
+        if _is_parameter(term.variable_1, p) ||
+           _is_parameter(term.variable_2, p)
+            continue
+        end
+        push!(H, (term.variable_1.value, term.variable_2.value))
+    end
+    return
+end
+
+function _append_sparse_hessian_structure!(
+    ::MOI.ScalarAffineFunction,
+    H,
+    ::Dict,
+)
+    return nothing
+end
+
+function _eval_sparse_hessian(
+    ∇²f::AbstractVector{T},
+    f::MOI.ScalarQuadraticFunction{T},
+    σ::T,
+    p::Dict{Int64,T},
+)::Int where {T}
+    i = 0
+    for term in f.quadratic_terms
+        if _is_parameter(term.variable_1, p) ||
+           _is_parameter(term.variable_2, p)
+            continue
+        end
+        i += 1
+        ∇²f[i] = term.coefficient * σ
+    end
+    return i
+end
+
+function _eval_sparse_hessian(
+    ∇²f::AbstractVector{T},
+    f::MOI.ScalarAffineFunction{T},
+    σ::T,
+    p::Dict{Int64,T},
+)::Int where {T}
+    return 0
+end
+
+Base.length(block::QPBlockData) = length(block.bound_type)
+
+function MOI.set(
+    block::QPBlockData{T},
+    ::MOI.ObjectiveFunction{F},
+    f::F,
+) where {T,F<:Union{MOI.VariableIndex,MOI.ScalarAffineFunction{T}}}
+    block.objective = convert(MOI.ScalarAffineFunction{T}, f)
+    block.objective_function_type = _function_info(f)
+    return
+end
+
+function MOI.set(
+    block::QPBlockData{T},
+    ::MOI.ObjectiveFunction{MOI.ScalarQuadraticFunction{T}},
+    f::MOI.ScalarQuadraticFunction{T},
+) where {T}
+    block.objective = f
+    block.objective_function_type = _function_info(f)
+    return
+end
+
+function MOI.get(block::QPBlockData{T}, ::MOI.ObjectiveFunctionType) where {T}
+    return _function_type_to_func(T, block.objective_function_type)
+end
+
+function MOI.get(block::QPBlockData{T}, ::MOI.ObjectiveFunction{F}) where {T,F}
+    return convert(F, block.objective)
+end
+
+function MOI.get(
+    block::QPBlockData{T},
+    ::MOI.ListOfConstraintTypesPresent,
+) where {T}
+    constraints = Set{Tuple{Type,Type}}()
+    for i in 1:length(block)
+        F = _function_type_to_func(T, block.function_type[i])
+        S = _bound_type_to_set(T, block.bound_type[i])
+        push!(constraints, (F, S))
+    end
+    return collect(constraints)
+end
+
+function MOI.is_valid(
+    block::QPBlockData{T},
+    ci::MOI.ConstraintIndex{F,S},
+) where {
+    T,
+    F<:Union{MOI.ScalarAffineFunction{T},MOI.ScalarQuadraticFunction{T}},
+    S<:Union{MOI.LessThan{T},MOI.GreaterThan{T},MOI.EqualTo{T},MOI.Interval{T}},
+}
+    return 1 <= ci.value <= length(block)
+end
+
+function MOI.get(
+    block::QPBlockData{T},
+    ::MOI.ListOfConstraintIndices{F,S},
+) where {
+    T,
+    F<:Union{MOI.ScalarAffineFunction{T},MOI.ScalarQuadraticFunction{T}},
+    S<:Union{MOI.LessThan{T},MOI.GreaterThan{T},MOI.EqualTo{T},MOI.Interval{T}},
+}
+    ret = MOI.ConstraintIndex{F,S}[]
+    for i in 1:length(block)
+        if _bound_type_to_set(T, block.bound_type[i]) != S
+            continue
+        elseif _function_type_to_func(T, block.function_type[i]) != F
+            continue
+        end
+        push!(ret, MOI.ConstraintIndex{F,S}(i))
+    end
+    return ret
+end
+
+function MOI.get(
+    block::QPBlockData{T},
+    ::MOI.NumberOfConstraints{F,S},
+) where {
+    T,
+    F<:Union{MOI.ScalarAffineFunction{T},MOI.ScalarQuadraticFunction{T}},
+    S<:Union{MOI.LessThan{T},MOI.GreaterThan{T},MOI.EqualTo{T},MOI.Interval{T}},
+}
+    return length(MOI.get(block, MOI.ListOfConstraintIndices{F,S}()))
+end
+
+function MOI.add_constraint(
+    block::QPBlockData{T},
+    f::Union{MOI.ScalarAffineFunction{T},MOI.ScalarQuadraticFunction{T}},
+    s::Union{MOI.LessThan{T},MOI.GreaterThan{T},MOI.EqualTo{T},MOI.Interval{T}},
+) where {T}
+    push!(block.constraints, f)
+    bound_type, l, u = _set_info(s)
+    push!(block.g_L, l)
+    push!(block.g_U, u)
+    push!(block.mult_g, nothing)
+    push!(block.bound_type, bound_type)
+    push!(block.function_type, _function_info(f))
+    return MOI.ConstraintIndex{typeof(f),typeof(s)}(length(block.bound_type))
+end
+
+function MOI.get(
+    block::QPBlockData{T},
+    ::MOI.ConstraintFunction,
+    c::MOI.ConstraintIndex{F,S},
+) where {T,F,S}
+    return convert(F, block.constraints[c.value])
+end
+
+function MOI.get(
+    block::QPBlockData{T},
+    ::MOI.ConstraintSet,
+    c::MOI.ConstraintIndex{F,S},
+) where {T,F,S}
+    row = c.value
+    if block.bound_type[row] == _kBoundTypeEqualTo
+        return MOI.EqualTo(block.g_L[row])::S
+    elseif block.bound_type[row] == _kBoundTypeLessThan
+        return MOI.LessThan(block.g_U[row])::S
+    elseif block.bound_type[row] == _kBoundTypeGreaterThan
+        return MOI.GreaterThan(block.g_L[row])::S
+    else
+        @assert block.bound_type[row] == _kBoundTypeInterval
+        return MOI.Interval(block.g_L[row], block.g_U[row])::S
+    end
+end
+
+function MOI.set(
+    block::QPBlockData{T},
+    ::MOI.ConstraintSet,
+    c::MOI.ConstraintIndex{F,MOI.LessThan{T}},
+    set::MOI.LessThan{T},
+) where {T,F}
+    block.g_U[c.value] = set.upper
+    return
+end
+
+function MOI.set(
+    block::QPBlockData{T},
+    ::MOI.ConstraintSet,
+    c::MOI.ConstraintIndex{F,MOI.GreaterThan{T}},
+    set::MOI.GreaterThan{T},
+) where {T,F}
+    block.g_L[c.value] = set.lower
+    return
+end
+
+function MOI.set(
+    block::QPBlockData{T},
+    ::MOI.ConstraintSet,
+    c::MOI.ConstraintIndex{F,MOI.EqualTo{T}},
+    set::MOI.EqualTo{T},
+) where {T,F}
+    block.g_L[c.value] = set.value
+    block.g_U[c.value] = set.value
+    return
+end
+
+function MOI.set(
+    block::QPBlockData{T},
+    ::MOI.ConstraintSet,
+    c::MOI.ConstraintIndex{F,MOI.Interval{T}},
+    set::MOI.Interval{T},
+) where {T,F}
+    block.g_L[c.value] = set.lower
+    block.g_U[c.value] = set.upper
+    return
+end
+
+function MOI.get(
+    block::QPBlockData{T},
+    ::MOI.ConstraintDualStart,
+    c::MOI.ConstraintIndex{F,S},
+) where {T,F,S}
+    return block.mult_g[c.value]
+end
+
+function MOI.set(
+    block::QPBlockData{T},
+    ::MOI.ConstraintDualStart,
+    c::MOI.ConstraintIndex{F,S},
+    value,
+) where {T,F,S}
+    block.mult_g[c.value] = value
+    return
+end
+
+function MOI.eval_objective(
+    block::QPBlockData{T},
+    x::AbstractVector{T},
+) where {T}
+    return _eval_function(block.objective, x, block.parameters)
+end
+
+function MOI.eval_objective_gradient(
+    block::QPBlockData{T},
+    ∇f::AbstractVector{T},
+    x::AbstractVector{T},
+) where {T}
+    ∇f .= zero(T)
+    _eval_dense_gradient(∇f, block.objective, x, block.parameters)
+    return
+end
+
+function MOI.eval_constraint(
+    block::QPBlockData{T},
+    g::AbstractVector{T},
+    x::AbstractVector{T},
+) where {T}
+    for (i, constraint) in enumerate(block.constraints)
+        g[i] = _eval_function(constraint, x, block.parameters)
+    end
+    return
+end
+
+function MOI.jacobian_structure(block::QPBlockData)
+    J = Tuple{Int,Int}[]
+    for (row, constraint) in enumerate(block.constraints)
+        _append_sparse_gradient_structure!(constraint, J, row, block.parameters)
+    end
+    return J
+end
+
+function MOI.eval_constraint_jacobian(
+    block::QPBlockData{T},
+    J::AbstractVector{T},
+    x::AbstractVector{T},
+) where {T}
+    i = 0
+    for constraint in block.constraints
+        ∇f = view(J, (i+1):length(J))
+        i += _eval_sparse_gradient(∇f, constraint, x, block.parameters)
+    end
+    return
+end
+
+function MOI.hessian_lagrangian_structure(block::QPBlockData)
+    H = Tuple{Int,Int}[]
+    _append_sparse_hessian_structure!(block.objective, H, block.parameters)
+    for constraint in block.constraints
+        _append_sparse_hessian_structure!(constraint, H, block.parameters)
+    end
+    return H
+end
+
+function MOI.eval_hessian_lagrangian(
+    block::QPBlockData{T},
+    H::AbstractVector{T},
+    x::AbstractVector{T},
+    σ::T,
+    μ::AbstractVector{T},
+) where {T}
+    i = _eval_sparse_hessian(H, block.objective, σ, block.parameters)
+    for (row, constraint) in enumerate(block.constraints)
+        ∇²f = view(H, (i+1):length(H))
+        i += _eval_sparse_hessian(∇²f, constraint, μ[row], block.parameters)
+    end
+    return
+end
+
+# The product functions below ACCUMULATE into their output vector, so that
+# the contributions of several blocks (for example, the QP block, the
+# vector-nonlinear-oracle constraints, and an `MOI.AbstractNLPEvaluator`) can
+# be composed into the same output. This is why they are not methods of the
+# corresponding `MOI.eval_...` functions, whose contract is to store the
+# result: `QPBlockData` is not an `MOI.AbstractNLPEvaluator`, so it does not
+# have to define the same interface as evaluators.
+
+function _add_Jv_product(
+    f::MOI.ScalarAffineFunction{T},
+    y::AbstractVector{T},
+    x::AbstractVector{T},
+    w::AbstractVector{T},
+    p::Dict{Int64,T},
+    i::Int,
+)::Nothing where {T}
+    for term in f.terms
+        if !_is_parameter(term.variable, p)
+            y[i] += term.coefficient * w[term.variable.value]
+        end
+    end
+    return
+end
+
+function _add_Jv_product(
+    f::MOI.ScalarQuadraticFunction{T},
+    y::AbstractVector{T},
+    x::AbstractVector{T},
+    w::AbstractVector{T},
+    p::Dict{Int64,T},
+    i::Int,
+)::Nothing where {T}
+    for term in f.affine_terms
+        if !_is_parameter(term.variable, p)
+            y[i] += term.coefficient * w[term.variable.value]
+        end
+    end
+    for term in f.quadratic_terms
+        if !_is_parameter(term.variable_1, p)
+            v = _value(term.variable_2, x, p)
+            y[i] += term.coefficient * v * w[term.variable_1.value]
+        end
+        if term.variable_1 != term.variable_2 &&
+           !_is_parameter(term.variable_2, p)
+            v = _value(term.variable_1, x, p)
+            y[i] += term.coefficient * v * w[term.variable_2.value]
+        end
+    end
+    return
+end
+
+function _add_Jtv_product(
+    f::MOI.ScalarAffineFunction{T},
+    y::AbstractVector{T},
+    x::AbstractVector{T},
+    w::AbstractVector{T},
+    p::Dict{Int64,T},
+    i::Int,
+)::Nothing where {T}
+    for term in f.terms
+        if !_is_parameter(term.variable, p)
+            y[term.variable.value] += term.coefficient * w[i]
+        end
+    end
+    return
+end
+
+function _add_Jtv_product(
+    f::MOI.ScalarQuadraticFunction{T},
+    y::AbstractVector{T},
+    x::AbstractVector{T},
+    w::AbstractVector{T},
+    p::Dict{Int64,T},
+    i::Int,
+)::Nothing where {T}
+    for term in f.affine_terms
+        if !_is_parameter(term.variable, p)
+            y[term.variable.value] += term.coefficient * w[i]
+        end
+    end
+    for term in f.quadratic_terms
+        if !_is_parameter(term.variable_1, p)
+            v = _value(term.variable_2, x, p)
+            y[term.variable_1.value] += term.coefficient * v * w[i]
+        end
+        if term.variable_1 != term.variable_2 &&
+           !_is_parameter(term.variable_2, p)
+            v = _value(term.variable_1, x, p)
+            y[term.variable_2.value] += term.coefficient * v * w[i]
+        end
+    end
+    return
+end
+
+function _add_Hv_product(
+    f::MOI.ScalarQuadraticFunction{T},
+    H::AbstractVector{T},
+    x::AbstractVector{T},
+    v::AbstractVector{T},
+    λ::T,
+    p::Dict{Int64,T},
+)::Nothing where {T}
+    for term in f.quadratic_terms
+        if _is_parameter(term.variable_1, p) ||
+           _is_parameter(term.variable_2, p)
+            continue
+        end
+        i, j = term.variable_1.value, term.variable_2.value
+        H[i] += λ * term.coefficient * v[j]
+        if i != j
+            H[j] += λ * term.coefficient * v[i]
+        end
+    end
+    return
+end
+
+function _add_Hv_product(
+    ::MOI.ScalarAffineFunction{T},
+    H::AbstractVector{T},
+    x::AbstractVector{T},
+    v::AbstractVector{T},
+    λ::T,
+    p::Dict{Int64,T},
+) where {T}
+    return nothing
+end
+
+# These are used to add the QP contribution on top of the NL contribution.
+
+"""
+    add_constraint_jacobian_product(
+        block::QPBlockData{T},
+        y::AbstractVector{T},
+        x::AbstractVector{T},
+        w::AbstractVector{T},
+    )::Nothing where {T}
+
+Add to `y` the product of the Jacobian of the constraints of `block` at `x`
+with `w`.
+
+Unlike [`MOI.eval_constraint_jacobian_product`](@ref), this function
+accumulates into `y` instead of storing the result, so that the contributions
+of several blocks can be composed: the caller is responsible for zeroing `y`
+before the first contribution.
+"""
+function add_constraint_jacobian_product(
+    block::QPBlockData{T},
+    y::AbstractVector{T},
+    x::AbstractVector{T},
+    w::AbstractVector{T},
+) where {T}
+    for (i, constraint) in enumerate(block.constraints)
+        _add_Jv_product(constraint, y, x, w, block.parameters, i)
+    end
+    return
+end
+
+"""
+    add_constraint_jacobian_transpose_product(
+        block::QPBlockData{T},
+        y::AbstractVector{T},
+        x::AbstractVector{T},
+        w::AbstractVector{T},
+    )::Nothing where {T}
+
+Add to `y` the product of the transpose of the Jacobian of the constraints of
+`block` at `x` with `w`.
+
+Unlike [`MOI.eval_constraint_jacobian_transpose_product`](@ref), this
+function accumulates into `y` instead of storing the result, so that the
+contributions of several blocks can be composed: the caller is responsible
+for zeroing `y` before the first contribution.
+"""
+function add_constraint_jacobian_transpose_product(
+    block::QPBlockData{T},
+    y::AbstractVector{T},
+    x::AbstractVector{T},
+    w::AbstractVector{T},
+) where {T}
+    for (i, constraint) in enumerate(block.constraints)
+        _add_Jtv_product(constraint, y, x, w, block.parameters, i)
+    end
+    return
+end
+
+"""
+    add_hessian_lagrangian_product(
+        block::QPBlockData{T},
+        H::AbstractVector{T},
+        x::AbstractVector{T},
+        v::AbstractVector{T},
+        σ::T,
+        μ::AbstractVector{T},
+    )::Nothing where {T}
+
+Add to `H` the product of the Hessian of the Lagrangian of `block` at `x`,
+with objective weight `σ` and constraint weights `μ`, with `v`.
+
+Unlike [`MOI.eval_hessian_lagrangian_product`](@ref), this function
+accumulates into `H` instead of storing the result, so that the contributions
+of several blocks can be composed: the caller is responsible for zeroing `H`
+before the first contribution.
+"""
+function add_hessian_lagrangian_product(
+    block::QPBlockData{T},
+    H::AbstractVector{T},
+    x::AbstractVector{T},
+    v::AbstractVector{T},
+    σ::T,
+    μ::AbstractVector{T},
+) where {T}
+    _add_Hv_product(block.objective, H, x, v, σ, block.parameters)
+    for (i, constraint) in enumerate(block.constraints)
+        _add_Hv_product(constraint, H, x, v, μ[i], block.parameters)
+    end
+    return
+end