Quick start

The primary purpose of ConstraintTrees.jl is to make the representation of constraint systems neat, and thus make their manipulation easy and high-level. In short, the package abstracts the users from keeping track of variable and constraint indexes in matrix form, and gives a nice data structure that describes the system, while keeping all variable allocation and constraint organization completely implicit.

Here we demonstrate the absolutely basic concepts on the "field allocation" problem.

The problem: Field area allocation

Suppose we have 100 square kilometers of field, 500 kilos of fertilizer and 300 kilos of insecticide. We also have a practically infinite supply of wheat and barley seeds. If we decide to sow barley, we can make 550🪙 per square kilometer of harvest; if we decide to go with wheat instead, we can make 350🪙. Unfortunately each square kilometer of wheat requires 6 kilos of fertilizer, and 1 kilo of insecticide, whereas each square kilometer of barley requires 2 kilos of fertilizer and 4 kilos of insecticide, because insects love barley.

How much of our fields should we allocate to wheat and barley to maximize our profit?

Field area allocation with ConstraintTrees.jl

Let's import the package and start constructing the problem:

import ConstraintTrees as C

Let's name our system s. We first need a few variables:

s = C.variables(keys = [:wheat, :barley])

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 2 elements:
  :barley => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#))
  :wheat  => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#))

With ConstraintTrees.jl, we can (and want to!) label everything very nicely – the constraint trees are essentially directory structures, so one can prefix everything with symbols to put it into nice directories, e.g. as such:

:area^s

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 1 element:
  :area => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)

To be absolutely realistic, we also want to make sure that all areas are non-negative. To demonstrate how to do that nicely from the start, we rather re-do the constraints with an appropriate interval bound:

s = :area^C.variables(keys = [:wheat, :barley], bounds = C.Between(0, Inf))

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 1 element:
  :area => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)

Constraint trees can be browsed using dot notation, or just like dictionaries:

s.area

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 2 elements:
  :barley => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#), ConstraintTrees.Between(0.0, Inf))
  :wheat  => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#), ConstraintTrees.Between(0.0, Inf))

s[:area].barley

ConstraintTrees.Constraint(ConstraintTrees.LinearValue([2], [1.0]), ConstraintTrees.Between(0.0, Inf))

(For convenience in some cases, string indexes are also supported:)

s["area"]["barley"]

ConstraintTrees.Constraint(ConstraintTrees.LinearValue([2], [1.0]), ConstraintTrees.Between(0.0, Inf))

Now let's start rewriting the problem into the constraint-tree-ish description. First, we only have 100 square kilometers of area:

total_area = s.area.wheat.value + s.area.barley.value

total_area_constraint = C.Constraint(total_area, (0, 100))

ConstraintTrees.Constraint(ConstraintTrees.LinearValue([1, 2], [1.0, 1.0]), ConstraintTrees.Between(0.0, 100.0))

We can add any kind of constraint into the existing constraint trees by "merging" multiple trees with operator *:

s *= :total_area^total_area_constraint

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 2 elements:
  :area       => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :total_area => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#), ConstraintTrees.Between(0.0, 100.0))

Now let's add constraints for resources. We can create whole ConstraintTree structures like dictionaries in place, as follows:

s *=
    :resources^C.ConstraintTree(
        :fertilizer =>
            C.Constraint(s.area.wheat.value * 6 + s.area.barley.value * 2, (0, 500)),
        :insecticide =>
            C.Constraint(s.area.wheat.value * 1 + s.area.barley.value * 4, (0, 300)),
    )

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 3 elements:
  :area       => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :resources  => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :total_area => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#), ConstraintTrees.Between(0.0, 100.0))

We can also represent the expected profit as a constraint (although we do not need to actually put a constraint bound there):

s *= :profit^C.Constraint(s.area.wheat.value * 350 + s.area.barley.value * 550)

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 4 elements:
  :area       => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :profit     => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#))
  :resources  => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :total_area => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#), ConstraintTrees.Between(0.0, 100.0))

Solving the system with JuMP

We can now take the structure of the constraint tree, translate it to any suitable linear optimizer interface, and have it solved. For popular reasons we choose JuMP with GLPK – the code is left uncommented here as-is; see the other examples for a slightly more detailed explanation:

import JuMP
function optimized_vars(cs::C.ConstraintTree, objective::C.LinearValue, optimizer)
    model = JuMP.Model(optimizer)
    JuMP.@variable(model, x[1:C.var_count(cs)])
    JuMP.@objective(model, JuMP.MAX_SENSE, C.substitute(objective, x))
    C.traverse(cs) do c
        b = c.bound
        if b isa C.EqualTo
            JuMP.@constraint(model, C.substitute(c.value, x) == b.equal_to)
        elseif b isa C.Between
            val = C.substitute(c.value, x)
            isinf(b.lower) || JuMP.@constraint(model, val >= b.lower)
            isinf(b.upper) || JuMP.@constraint(model, val <= b.upper)
        end
    end
    JuMP.optimize!(model)
    JuMP.value.(model[:x])
end

import GLPK
optimal_variable_assignment = optimized_vars(s, s.profit.value, GLPK.Optimizer)

2-element Vector{Float64}:
 33.333333333333336
 66.66666666666667

This gives us the optimized variable values! If we cared to remember what they stand for, we might already know how much barley to sow. On the other hand, the main point of ConstraintTree.jl is that one should not be forced to remember things like variable ordering and indexes, or be forced to manually calculate how much money we actually make or how much fertilizer is going to be left – instead, we can simply feed the variable values back to the constraint tree, and get a really good overview of all values in our constrained system:

optimal_s = C.substitute_values(s, optimal_variable_assignment)

ConstraintTrees.Tree{Float64} with 4 elements:
  :area       => ConstraintTrees.Tree{Float64}(#= 2 elements =#)
  :profit     => 48333.3
  :resources  => ConstraintTrees.Tree{Float64}(#= 2 elements =#)
  :total_area => 100.0

Browsing the result

optimal_s is now like the original constraint tree, just the contents are "plain old values" instead of the constraints as above. Thus we can easily see our profit:

optimal_s.profit

48333.33333333334

The occupied area for each crop:

optimal_s.area

ConstraintTrees.Tree{Float64} with 2 elements:
  :barley => 66.6667
  :wheat  => 33.3333

The consumed resources:

optimal_s.resources

ConstraintTrees.Tree{Float64} with 2 elements:
  :fertilizer  => 333.333
  :insecticide => 300.0

Increasing the complexity

A crucial property of constraint trees is that the users do not need to care about what kind of value they are manipulating – no matter if something is a variable or a derived value, the code that works with it is the same. For example, we can use the actual prices for our resources (30🪙 and 110🪙 for a kilo of fertilizer and insecticide, respectively) to make a corrected profit:

s *=
    :actual_profit^C.Constraint(
        s.profit.value - 30 * s.resources.fertilizer.value -
        110 * s.resources.insecticide.value,
    )

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 5 elements:
  :actual_profit => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#))
  :area          => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :profit        => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#))
  :resources     => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :total_area    => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#), ConstraintTrees.Between(0.0, 100.0))

Is the result going to change if we optimize for the corrected profit?

realistically_optimal_s =
    C.substitute_values(s, optimized_vars(s, s.actual_profit.value, GLPK.Optimizer))

ConstraintTrees.Tree{Float64} with 5 elements:
  :actual_profit => 5750.0
  :area          => ConstraintTrees.Tree{Float64}(#= 2 elements =#)
  :profit        => 40000.0
  :resources     => ConstraintTrees.Tree{Float64}(#= 2 elements =#)
  :total_area    => 100.0

realistically_optimal_s.area

ConstraintTrees.Tree{Float64} with 2 elements:
  :barley => 25.0
  :wheat  => 75.0

Combining constraint systems: Let's have a factory!

The second crucial property of constraint trees is the ability to easily combine different constraint systems into one. Let's pretend we also somehow obtained a food factory that produces malty sweet bread and wheaty weizen-style beer, with various extra consumptions of water and heat for each of the products. For simplicity, let's just create the corresponding constraint system (f as a factory) here:

f = :products^C.variables(keys = [:bread, :weizen], bounds = C.Between(0, Inf))
f *= :profit^C.Constraint(25 * f.products.weizen.value + 35 * f.products.bread.value)

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 2 elements:
  :products => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :profit   => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#))

We can make the constraint systems more complex by adding additional variables. To make sure the variables do not "conflict", one must use the + operator. While constraint systems combined with * always share variables, constraint systems combined with + are independent.

f += :materials^C.variables(keys = [:wheat, :barley], bounds = C.Between(0, Inf))

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 3 elements:
  :materials => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :products  => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :profit    => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#))

How much resources are consumed by each product, with a limit on each:

f *=
    :resources^C.ConstraintTree(
        :heat => C.Constraint(
            5 * f.products.bread.value + 3 * f.products.weizen.value,
            (0, 1000),
        ),
        :water => C.Constraint(
            2 * f.products.bread.value + 10 * f.products.weizen.value,
            (0, 3000),
        ),
    )

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 4 elements:
  :materials => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :products  => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :profit    => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#))
  :resources => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)

How much raw materials are required for each product:

f *=
    :material_allocation^C.ConstraintTree(
        :wheat => C.Constraint(
            8 * f.products.bread.value + 2 * f.products.weizen.value -
            f.materials.wheat.value,
            0,
        ),
        :barley => C.Constraint(
            0.5 * f.products.bread.value + 10 * f.products.weizen.value -
            f.materials.barley.value,
            0,
        ),
    )

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 5 elements:
  :material_allocation => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :materials           => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :products            => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)
  :profit              => ConstraintTrees.Constraint(ConstraintTrees.LinearValue(#= ... =#))
  :resources           => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)

Having the two systems at hand, we can connect the factory "system" f to the field "system" s, making a compound system c as such:

c = :factory^f + :fields^s

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 2 elements:
  :factory => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 5 elements =#)
  :fields  => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 5 elements =#)

Operators for combining constraint trees

Always remember to use + instead of * when combining independent constraint systems. If we use *, the variables in both systems will become implicitly shared, which is rarely what one wants in the first place. Use * only if adding additional constraints to an existing system. As a rule of thumb, one can remember the boolean interpretation of * as "and" and of + as "or".

On a side note, the operator ^ was chosen mainly to match the algebraic view of the tree combination, and nicely fit into Julia's operator priority structure.

To actually connect the systems (which now exist as completely independent parts of s), let's add a transport – the barley and wheat produced on the fields is going to be the only barley and wheat consumed by the factory, thus their production and consumption must sum to net zero:

c *= :transport^C.zip(c.fields.area, c.factory.materials) do area, material
    C.Constraint(area.value - material.value, 0)
end

ConstraintTrees.Tree{ConstraintTrees.Constraint} with 3 elements:
  :factory   => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 5 elements =#)
  :fields    => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 5 elements =#)
  :transport => ConstraintTrees.Tree{ConstraintTrees.Constraint}(#= 2 elements =#)

High-level constraint tree manipulation

There is also a dedicated example with many more useful functions like zip above.

Finally, let's see how much money can we make from having the factory supported by our fields in total!

optimal_c =
    C.substitute_values(c, optimized_vars(c, c.factory.profit.value, GLPK.Optimizer))

ConstraintTrees.Tree{Float64} with 3 elements:
  :factory   => ConstraintTrees.Tree{Float64}(#= 5 elements =#)
  :fields    => ConstraintTrees.Tree{Float64}(#= 5 elements =#)
  :transport => ConstraintTrees.Tree{Float64}(#= 2 elements =#)

How much field area did we allocate?

optimal_c.fields.area

ConstraintTrees.Tree{Float64} with 2 elements:
  :barley => 25.0
  :wheat  => 75.0

How much of each of the products does the factory make in the end?

optimal_c.factory.products

ConstraintTrees.Tree{Float64} with 2 elements:
  :bread  => 8.86076
  :weizen => 2.05696

How much extra resources is consumed by the factory?

optimal_c.factory.resources

ConstraintTrees.Tree{Float64} with 2 elements:
  :heat  => 50.4747
  :water => 38.2911

And what is the factory profit in the end?

optimal_c.factory.profit

361.5506329113926

This page was generated using Literate.jl.