Complex number support

This page explains the complex-valued variables and constraints that JuMP supports. For a worked-example using these features, read the Example: quantum state discrimination tutorial.

Complex-valued variables

Create a complex-valued variable using ComplexPlane:

julia> model = Model();julia> @variable(model, x in ComplexPlane())real(x) + imag(x) im

Note that x is not a VariableRef; instead, it is an affine expression with Complex{Float64}-valued coefficients:

julia> typeof(x)GenericAffExpr{ComplexF64, VariableRef}

Behind the scenes, JuMP has created two real-valued variables, with names "real(x)" and "imag(x)":

julia> all_variables(model)2-element Vector{VariableRef}: real(x) imag(x)julia> name.(all_variables(model))2-element Vector{String}: "real(x)" "imag(x)"

Use the real and imag functions on x to return a real-valued affine expression representing each variable:

julia> typeof(real(x))AffExpr (alias for GenericAffExpr{Float64, GenericVariableRef{Float64}})julia> typeof(imag(x))AffExpr (alias for GenericAffExpr{Float64, GenericVariableRef{Float64}})

To create an anonymous variable, use the set keyword argument:

julia> model = Model();julia> x = @variable(model, set = ComplexPlane())_[1] + _[2] im

Complex-valued variable and start values bounds

Because complex-valued variables lack a total ordering, the definition of a variable bound for a complex-valued variable is ambiguous. If you pass a real- or complex-valued argument to keywords such as lower_bound, upper_bound, and start_value, JuMP will apply the real and imaginary parts to the associated real-valued variables.

julia> model = Model();julia> @variable(           model,           x in ComplexPlane(),           lower_bound = 1.0,           upper_bound = 2.0 + 3.0im,           start = 4im,       )real(x) + imag(x) imjulia> vars = all_variables(model)2-element Vector{VariableRef}: real(x) imag(x)julia> lower_bound.(vars)2-element Vector{Float64}: 1.0 0.0julia> upper_bound.(vars)2-element Vector{Float64}: 2.0 3.0julia> start_value.(vars)2-element Vector{Float64}: 0.0 4.0

You can modify the bounds and start values by passing imag(x) or real(x) to the appropriate function:

julia> set_lower_bound(imag(x), 2)julia> lower_bound(imag(x))2.0julia> delete_upper_bound(real(x))julia> has_upper_bound(real(x))falsejulia> set_start_value(imag(x), 3)julia> start_value(imag(x))3.0

Complex-valued equality constraints

JuMP reformulates complex-valued equality constraints into two real-valued constraints: one representing the real part, and one representing the imaginary part. Thus, complex-valued equality constraints can be solved any solver that supports the real-valued constraint type.

For example:

julia> model = Model(HiGHS.Optimizer);julia> set_silent(model)julia> @variable(model, x[1:2]);julia> @constraint(model, (1 + 2im) * x[1] + 3 * x[2] == 4 + 5im)(1 + 2im) x[1] + 3 x[2] = (4 + 5im)julia> optimize!(model)julia> value.(x)2-element Vector{Float64}: 2.5 0.5

is equivalent to

julia> model = Model(HiGHS.Optimizer);julia> set_silent(model)julia> @variable(model, x[1:2]);julia> @constraint(model, 1 * x[1] + 3 * x[2] == 4)  # real componentx[1] + 3 x[2] = 4julia> @constraint(model, 2 * x[1] == 5)  # imag component2 x[1] = 5julia> optimize!(model)julia> value.(x)2-element Vector{Float64}: 2.5 0.5

This also applies if the variables are complex-valued:

julia> model = Model(HiGHS.Optimizer);julia> set_silent(model)julia> @variable(model, x in ComplexPlane());julia> @constraint(model, (1 + 2im) * x + 3 * x == 4 + 5im)(4 + 2im) real(x) + (-2 + 4im) imag(x) = (4 + 5im)julia> optimize!(model)julia> value(x)1.3 + 0.6000000000000001im

which is equivalent to

julia> model = Model(HiGHS.Optimizer);julia> set_silent(model)julia> @variable(model, x_real);julia> @variable(model, x_imag);julia> @constraint(model, x_real - 2 * x_imag + 3 * x_real == 4)4 x_real - 2 x_imag = 4julia> @constraint(model, x_imag + 2 * x_real + 3 * x_imag == 5)2 x_real + 4 x_imag = 5julia> optimize!(model)julia> value(x_real) + value(x_imag) * im1.3 + 0.6000000000000001im

Hermitian PSD Cones

JuMP supports creating matrices where are Hermitian.

julia> model = Model();julia> @variable(model, H[1:3, 1:3] in HermitianPSDCone())3×3 LinearAlgebra.Hermitian{GenericAffExpr{ComplexF64, VariableRef}, Matrix{GenericAffExpr{ComplexF64, VariableRef}}}: real(H[1,1])                    …  real(H[1,3]) + imag(H[1,3]) im real(H[1,2]) - imag(H[1,2]) im     real(H[2,3]) + imag(H[2,3]) im real(H[1,3]) - imag(H[1,3]) im     real(H[3,3])

Behind the scenes, JuMP has created nine real-valued decision variables:

julia> all_variables(model)9-element Vector{VariableRef}: real(H[1,1]) real(H[1,2]) real(H[2,2]) real(H[1,3]) real(H[2,3]) real(H[3,3]) imag(H[1,2]) imag(H[1,3]) imag(H[2,3])

and a Vector{VariableRef}-in-MOI.HermitianPositiveSemidefiniteConeTriangle constraint:

julia> num_constraints(model, Vector{VariableRef}, MOI.HermitianPositiveSemidefiniteConeTriangle)1

The MOI.HermitianPositiveSemidefiniteConeTriangle set can be efficiently bridged to MOI.PositiveSemidefiniteConeTriangle, so it can be solved by any solver that supports PSD constraints.

Each element of H is an affine expression with Complex{Float64}-valued coefficients:

julia> typeof(H[1, 1])GenericAffExpr{ComplexF64, VariableRef}julia> typeof(H[2, 1])GenericAffExpr{ComplexF64, VariableRef}

Start values

When setting the start value, you must be careful to set only the upper triangle of real variables, and the upper triangle excluding the diagonal of imaginary variables:

julia> import LinearAlgebrajulia> function set_hermitian_start(           H::LinearAlgebra.Hermitian,           start::LinearAlgebra.Hermitian,       )           for j in 1:size(H, 2), i in 1:j               set_start_value(real(H[i, j]), real(start[i, j]))               if i < j                   set_start_value(imag(H[i, j]), imag(start[i, j]))               end           end           return       endset_hermitian_start (generic function with 1 method)julia> H0 = LinearAlgebra.Hermitian(           [1 (2+3im) (5+6im); (2-3im) 4 (7+8im); (5-6im) (7-8im) 9],       )3×3 LinearAlgebra.Hermitian{Complex{Int64}, Matrix{Complex{Int64}}}: 1+0im  2+3im  5+6im 2-3im  4+0im  7+8im 5-6im  7-8im  9+0imjulia> set_hermitian_start(H, H0)julia> start_value.(all_variables(model))9-element Vector{Float64}: 1.0 2.0 4.0 5.0 7.0 9.0 3.0 6.0 8.0

Hermitian PSD constraints

The HermitianPSDCone can also be used in the @constraint macro:

julia> model = Model();julia> @variable(model, x[1:2])2-element Vector{VariableRef}: x[1] x[2]julia> import LinearAlgebrajulia> H = LinearAlgebra.Hermitian([x[1] 1im; -1im -x[2]])2×2 LinearAlgebra.Hermitian{GenericAffExpr{ComplexF64, VariableRef}, Matrix{GenericAffExpr{ComplexF64, VariableRef}}}: x[1]  im -im   -x[2]julia> @constraint(model, H in HermitianPSDCone())[x[1]  im -im   -x[2]] ∈ HermitianPSDCone()