Calcoli scientifici con Julia/Massimizzazione profitto

La massimizzazione del profitto è un classico esempio di programmazione lineare che nella sua forma standard si esprime in 3 parti:

• Una funzione lineare da massimizzare:

$f(x_{1},x_{2})=c_{1}x_{1}+c_{2}x_{2}$

• Problema dei vincoli della seguente forma:

${\begin{matrix}a_{11}x_{1}+a_{12}x_{2}&\leq b_{1}\\a_{21}x_{1}+a_{22}x_{2}&\leq b_{2}\\a_{31}x_{1}+a_{32}x_{2}&\leq b_{3}\\\end{matrix}}$

• Variabili non negative

${\begin{matrix}x_{1}\geq 0\\x_{2}\geq 0\end{matrix}}$

Il problema è di solito espresso in forma matriciale nel seguente modo:

$\max\{\,\mathbf {c} ^{\mathsf {T}}\mathbf {x} \mid \mathbf {x} \in \mathbb {R} ^{n}\land A\mathbf {x} \leq \mathbf {b} \land \mathbf {x} \geq 0\,\}$

Esempio di calcolo con Julia

Un'azienda produce tavoli e mobiletti che intende vendere al prezzo rispettivamente di 400€ e 200€. Per produrli occorrono 6 unità di legno, 2 ore/lavoratore e 100 minuti di tempo-macchina per i tavoli e 10 unità di legno, 4 ore/lavoratore e 40 minuti di tempo-macchina per i mobiletti. A disposizione l'azienda possiede giornalmente 4000 unità di legno, 250 lavoratori per 8 ore al giorno cioè 2000 ore/lavoratore e 50 macchine che funzionano per 16 ore giornaliere cioè 48000 ore/macchina.

Innanzitutto occorre installare i pacchetti JuMP ed HiGHS.

Da REPL digitare:

using Pkg
Pkg.add("HiGHS")

using Pkg
Pkg.add("JuMP")

Da Jupyter in due celle di codice digitare :

] add HiGHS

] add JuMP

e poi cliccare su Run. Poi eseguire il codice per creare il modello di programmazione lineare:

using JuMP
using HiGHS
model = Model(HiGHS.Optimizer)
@variable(model, x >= 0)
@variable(model, y >= 0)
@objective(model, Max, 200y + 400x)
@constraint(model, c1, 10y + 6x <= 4000)
@constraint(model, c2, 4y + 2x <= 2000)
@constraint(model, c3, 40y + 100x <= 48000)
print(model)
optimize!(model)

max Subject to 200𝑦+400𝑥
6𝑥+10𝑦≤4000
2𝑥+4𝑦≤2000
100𝑥+40𝑦≤48000
𝑥≥0
𝑦≥0

Running HiGHS 1.5.3 [date: 1970-01-01, git hash: 45a127b78]
Copyright (c) 2023 HiGHS under MIT licence terms
Presolving model
3 rows, 2 cols, 6 nonzeros
3 rows, 2 cols, 6 nonzeros
Presolve : Reductions: rows 3(-0); columns 2(-0); elements 6(-0) - Not reduced
Problem not reduced by presolve: solving the LP
Using EKK dual simplex solver - serial
 Iteration        Objective     Infeasibilities num(sum)
         0    -3.7499962007e+01 Ph1: 3(5.6875); Du: 2(37.5) 0s
         2     1.9789473684e+05 Pr: 0(0) 0s
Model   status      : Optimal
Simplex   iterations: 2
Objective value     :  1.9789473684e+05
HiGHS run time      :          0.03

println("termination_status = ",termination_status(model))
println("primal_status = ",primal_status(model))
println("dual_status = ",dual_status(model))
println("objective_value = ",objective_value(model))
println("x = ",value(x))
println("y = ",value(y))
println("c1 = ",shadow_price(c1))
println("c2 = ",shadow_price(c2))
println("c3 = ",shadow_price(c3))

termination_status = OPTIMAL
primal_status = FEASIBLE_POINT
dual_status = FEASIBLE_POINT
objective_value = 197894.73684210528
x = 421.05263157894734
y = 147.36842105263167
c1 = 5.2631578947368425
c2 = -0.0
c3 = 3.6842105263157894

cioè l'azienda ottiene un profitto massimo giornaliero di 197.894€ producendo 147 mobiletti e 421 tavoli.