Fix for constraint with unit coefficients

[blegat]

This solves an issue that @trulsf is facing in https://github.com/trulsf/UnitJuMP.jl/tree/tf/units_in_affexpr

Before the PR:

ERROR: DimensionError: 0.0 and 1.0 m are not dimensionally compatible.
Stacktrace:
 [1] +(x::Quantity{Float64, NoDims, Unitful.FreeUnits{…}}, y::Quantity{Float64, 𝐋, Unitful.FreeUnits{…}})
   @ Unitful ~/.julia/packages/Unitful/GYzMo/src/quantities.jl:137
 [2] +(x::Bool, y::Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}})
   @ Base ./promotion.jl:422
 [3] shift_constant(set::MathOptInterface.EqualTo{Bool}, offset::Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}})
   @ MathOptInterface.Utilities ~/.julia/packages/MathOptInterface/aJZbq/src/Utilities/sets.jl:80
 [4] build_constraint(::Function, expr::GenericAffExpr{Quantity{…}, VariableRef}, set::MathOptInterface.EqualTo{Bool})
   @ JuMP ~/.julia/dev/JuMP/src/macros/@constraint.jl:908
 [5] build_constraint(::Function, ::GenericAffExpr{Quantity{…}, VariableRef}, ::Zeros; kwargs::@Kwargs{})
   @ JuMP ~/.julia/dev/JuMP/src/macros/@constraint.jl:828
 [6] build_constraint(::Function, ::GenericAffExpr{Quantity{Float64, 𝐋, Unitful.FreeUnits{…}}, VariableRef}, ::Zeros)
   @ JuMP ~/.julia/dev/JuMP/src/macros/@constraint.jl:827
 [7] macro expansion
   @ ~/.julia/dev/JuMP/src/macros/@constraint.jl:172 [inlined]
 [8] macro expansion
   @ ~/.julia/dev/JuMP/src/macros.jl:393 [inlined]
 [9] top-level scope
   @ REPL[17]:1

After using MOI.EqualTo(MA.Zero()):

julia> @constraint(model, a == 1.0u"m")
ERROR: DimensionError: 1.0 m and 0.0 are not dimensionally compatible.
Stacktrace:
 [1] +(x::Quantity{Float64, 𝐋, Unitful.FreeUnits{…}}, y::Quantity{Float64, NoDims, Unitful.FreeUnits{…}})
   @ Unitful ~/.julia/packages/Unitful/GYzMo/src/quantities.jl:137
 [2] +(x::Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}, y::Float64)
   @ Base ./promotion.jl:422
 [3] shift_constant(set::MathOptInterface.EqualTo{Quantity{Float64, 𝐋, Unitful.FreeUnits{…}}}, offset::Float64)
   @ MathOptInterface.Utilities ~/.julia/packages/MathOptInterface/aJZbq/src/Utilities/sets.jl:80
 [4] model_convert(model::Model, set::MathOptInterface.EqualTo{Quantity{Float64, 𝐋, Unitful.FreeUnits{(m,), 𝐋, nothing}}})
   @ JuMP ~/.julia/dev/JuMP/src/macros.jl:306
 [5] model_convert(model::Model, con::ScalarConstraint{GenericAffExpr{…}, MathOptInterface.EqualTo{…}})
   @ JuMP ~/.julia/dev/JuMP/src/macros.jl:327
 [6] macro expansion
   @ ~/.julia/dev/JuMP/src/macros/@constraint.jl:172 [inlined]
 [7] macro expansion
   @ ~/.julia/dev/JuMP/src/macros.jl:393 [inlined]
 [8] top-level scope
   @ REPL[9]:1
Some type information was truncated. Use `show(err)` to see complete types.

After commenting the lines in model_convert

julia> a = GenericAffExpr{typeof(1.0u"m"),typeof(x)}(0.0u"m", x => 1.0u"m")
x

julia> @constraint(model, a == 1.0u"m")
x = 1.0 m

Requires jump-dev/MutableArithmetics.jl#294

[blegat]

Actually, the following requires no modification in JuMP nor MA:

using JuMP, Unitful
struct UnitShape{U} <: JuMP.AbstractShape
    unit::U
end

function JuMP._complex_convert(::Type{T}, u::Unitful.Quantity{U,A,B}) where {T,U,A,B}
    return convert(Unitful.Quantity{T,A,B}, u)
end

JuMP.reshape_vector(v, s::UnitShape) = v[] * s.unit
JuMP.reshape_set(::MOI.Zeros, shape::UnitShape) = MOI.EqualTo(zero(shape.unit))

function JuMP.build_constraint(
    ::Function,
    expr::GenericAffExpr{U},
    set::MOI.EqualTo,
) where {U<:Unitful.Quantity}
    return VectorConstraint([map_coefficients(c -> c.val, expr) - MOI.constant(set)], MOI.Zeros(1), UnitShape(oneunit(U)))
end

using JuMP, Unitful, HiGHS
model = Model(HiGHS.Optimizer)
@variable(model, x)
a = GenericAffExpr{typeof(1.0u"m"),typeof(x)}(0.0u"m", x => 1.0u"m")
@constraint(model, a == 1u"m")
optimize!(model)
value(x)

We have to use VectorConstraint because ScalarConstraint does not support shapes.

[trulsf]

The approach seems to work ok for constraints, but I can not see how to extend it to the objective function(s). I am not familiar with the use of shapes, but it seems to be designed for a different use and the code becomes a bit hard to understand. I think using bridges would be a more natural approach. It also opens up for the possibility of having (future) solvers that can accept units in their input. That said, I have so far been unable to get a bridge up and running for testing.

[trulsf]

I now have a working version with a constraint bridge (LessThan) and an objective bridge that strips off units at https://github.com/trulsf/UnitJuMP.jl/tree/tf/units_in_affexpr

Still a bit simplistic and primitive, but at least it seems doable to use bridges.
It requires a minor modification in MathOptInterface to allow numbers in LessThan:

struct LessThan{T<:Number} <: AbstractScalarSet

[trulsf]

To follow up a bit on the last comment. I see that a GenericModel has a type argument T<:Real. There are also similar requirements for e.g. LessThan and GreaterThan. Now Quantity<:Number, and it is therefore not possible to make a GenericModel{Quantity} which I guess would be the correct thing to do when using units. Is there a specific reason why coefficients are restricted to being Real?

Working with a GenericModel{Float64} works for the model transformation, but e.g. implementing the objective_value of the bridge will not allow the return of a Quantity.

[blegat]

If you have constraints with different units then what you be the type U in GenericModel{U} ? The type T in GenericModel{T} is the type of variables. Since a variable with units is constructed as a GenericAffExpr where the variable is still Float64 then it should be a GenericModel{Float64}.
The fact that we are subtyping Real is because we expect it to be something comparable.

The approach seems to work ok for constraints, but I can not see how to extend it to the objective function(s).

Currently, only VectorConstraints support shapes but this is mostly because we only had a use case for this at the moment. We could add shapes for scalar constraints as well as objective now that we have a use case for it.

I am not familiar with the use of shapes, but it seems to be designed for a different use and the code becomes a bit hard to understand.

I think this is the use case it is designed for. The MOI backend of a GenericModel{T} expects variables to have values of type T, scalar functions to evaluate to values of type T and vector functions to evaluate to values of type Vector{T}. We sometimes want the user to work with other types. For instance, we want to support matrices for semidefinite programming and polynomials in PolyJuMP. The conversion from these types to a Vector{T} is done with shapes and JuMP.vectorize and JuMP.reshape_vector. For scalar functions and objective, vectorize and reshape_vector are probably not the right terminology but we could extended it anyway to the scalar case.

I think using bridges would be a more natural approach. It also opens up for the possibility of having (future) solvers that can accept units in their input.

We usually want to create new MOI sets in order for sets to support them natively. But for units, what would be the advantage of working with units ? I don't see how that would speed up the solving process.