# Performance of Linear Programming Problem

**URL:** https://discourse.julialang.org/t/performance-of-linear-programming-problem/74569
**Category:** Optimization (Mathematical)
**Tags:** jump, performance, gurobi, linear-programming
**Created:** [January 13, 2022, 6:17pm UTC](https://discourse.julialang.org/t/performance-of-linear-programming-problem/74569 "2022-01-13T18:17:06Z")
**Posts on this page:** 2
**Page:** 1

<div class="post-metadata">

### Author: ![jmcastro2109](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jmcastro2109/32/38427_2.png) [@jmcastro2109](https://discourse.julialang.org/u/jmcastro2109)
#### Post date: [January 13, 2022, 6:17pm UTC](https://discourse.julialang.org/t/performance-of-linear-programming-problem/74569/1 "2022-01-13T18:17:06Z")

</div>

Hi everybody! I am trying to improve the performance of a conventional Linear Programming Problem which I am solving by calling Gurobi through JuMP. A MWE is the following

```julia-auto
using JuMP, Ipopt, Optim, LinearAlgebra, Random, Distributions, Gurobi, BenchmarkTools

Ι = 500
J = 50

const GUROBI_ENV = Gurobi.Env()

function solve_model()
	model = Model(with_optimizer(Gurobi.Optimizer, GUROBI_ENV); add_bridges = false)
end

#Create a vector of productivities

z=exp.(rand(Normal(0,1),Ι))

# Create a Matrix of Distances

function distmat(J,Ι)
    Distances = zeros(J,Ι)
    Random.seed!(7777)
    coordinates_market = 100*rand(J,2)
    coordinates_plant = 100*rand(Ι,2)
    for j = 1:J, l=1:Ι
        Distances[j,l] = sqrt((coordinates_market[j,1]-coordinates_plant[l,1])^2+(coordinates_market[j,2]-coordinates_plant[l,2])^2)
    end
    return 1 .+ Distances./100
end

τ = distmat(J,Ι)

function solving_min_constraints_primal(Cap,Q,τ,z)
    (J,Ι) = size(τ)
    primal_capacity = solve_model()
    set_silent(primal_capacity)
    #set_optimizer_attribute(primal_capacity, "nlp_scaling_method", "none")
    @variable(primal_capacity, x[1:J,1:Ι] >= 0)
    @objective(primal_capacity, Min, sum((τ[:,i]'*x[:,i])/z[i] for i=1:Ι) )

    for j=1:J
		@constraint(primal_capacity, sum(x[j,:]) == Q[j])
	end

    for i=1:Ι
		@constraint(primal_capacity, τ[:,i]'*x[:,i] <= Cap[i])
	end

    optimize!(primal_capacity)

    
    return objective_value(primal_capacity), value.(x)

end 

res_primal = @btime solving_min_constraints_primal(ones(Ι).*1.1,ones(J),τ,z)

92.057 ms (701849 allocations: 56.02 MiB)

```

Is there any performance recommendation to improve the speed of this code? Is there a way to reduce the number of allocations as much as possible? Are there any gains from vectorizing stuff and writing things in terms of matrices (my own attempts seem to suggest that no)? Although, I have not explored if there is an improvement related to the Sparsity of the constraints.

Any help would be very much appreciated!

---

<div class="post-metadata">

### Author: ![miles.lubin](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/miles.lubin/32/279_2.png) [@miles.lubin](https://discourse.julialang.org/u/miles.lubin)
#### Post date: [January 13, 2022, 8:03pm UTC](https://discourse.julialang.org/t/performance-of-linear-programming-problem/74569/2 "2022-01-13T20:03:10Z")

</div>

> Are there any gains from vectorizing stuff and writing things in terms of matrices (my own attempts seem to suggest that no)?

No

Two suggestions:

1. Use profiling to identify where the time is spent. You can use `JuMP.solve_time` to query the time that Gurobi reports spending on actually solving the problem. If this number is a majority of the 92ms, then there’s no reason to try optimizing JuMP more.

2. Try [`JuMP.direct_model`](https://jump.dev/JuMP.jl/dev/reference/models/#JuMP.direct_model) instead of `JuMP.Model`. This removes a translation layer at the expense of less flexibility in modeling, which the above code doesn’t need.