::install("https://ampl.com/dl/API/rAMPL.tar.gz", repos=NULL, INSTALL_opts=c("--no-multiarch", "--no-staged-install")) renv
Chapter 1 - AMPL Book
Hi Folks!
Since my bachelor degree I always wanted to resolved the exercises of AMPL’s book, In my course I studied with AMPL and in my master I also used AMPL, I find it very easy to get into operation research world!.
Have fun with me during solving some exercises from the book, I have used the API’s that AMPL is giving us, so some problems are solved using python and others with R.
I tried to keep the original model and data in different files and folders. When a model requires to be changed, I do it using the ampl API.
0. Installing rAMPL and amplpy
rAMPL
amplpy
In order to use the free license of ampl community, it is required to input the license number after –uuid.
```{r}
reticulate::py_install("jupyter")
reticulate::py_install("amplpy")
reticulate::py_install("pandas")
system2(reticulate::py_exe(), c("-m", "amplpy.modules", "install", "highs"))
system2(reticulate::py_exe(), c("-m", "amplpy.modules", "run", "amplkey", "activate", "--UUID"))
system2(reticulate::py_exe(), c("-m", "amplpy.modules", "run", "ampl", "-vvq"))
```
1. Optimize advertising campaigns (rAMPL)
The objective of this exercise is to chose the number of units to purchase in each channel to advertise. Let’s figure out together how to do it.
Setting a new environment for make the conection through the API with AMPL
a. Subjetct to budget
# creating vectors with information
<- c("tv_mins", "magazine_pages")
units <- c(20000, 10000)
cost <- c(1.8, 1)
reach <- c(10, 1)
min_units
# Loading class and setting solver
<- new(AMPL, env)
ampl $setOption("solver","HiGHS")
ampl
$eval("param p_budget := 1000000;")
ampl# Loading model
$read("models/1.1 Advertising_campaigns.mod")
ampl
# Setting Data
$setData(
ampldata.frame(
units = units,
cost = cost,
reach = reach,
min_units = min_units
1, "units"
),
) # Formulation
$eval("var buy{u in units} >= min_units[u];")
ampl$eval("maximize Audience: sum {u in units} buy[u] * reach[u];")
ampl$eval("subject to budget: sum {u in units} buy[u] * cost[u] <= p_budget;")
ampl
$solve() ampl
HiGHS 1.7.0: HiGHS 1.7.0: optimal solution; objective 98
0 simplex iterations
0 barrier iterations
print(ampl$getVariable("buy")$getValues())
index0 buy.val
1 magazine_pages 80
2 tv_mins 10
$close() ampl
b. Adding constraint time needed to build campaign
<- c("tv_mins", "magazine_pages")
units <- c(20000, 10000)
cost <- c(1.8, 1)
reach <- c(10, 1)
min_units <- c(1, 3)
person_weeks
# loading class and setting solver
<- new(AMPL, env)
ampl $setOption("solver","highs")
ampl
# Adding parameters
$eval("param p_budget := 1000000;")
ampl$eval("param max_person_weeks := 100;")
ampl# loading model
$read("models/1.1 advertising_campaigns.mod") # read model located in folder models
ampl
# setting data
$setData(
ampldata.frame(
units = units,
cost = cost,
reach = reach,
min_units = min_units,
person_weeks = person_weeks
1, "units"
),
) # formulation
$eval("var buy{u in units} >= min_units[u];")
ampl$eval("maximize audience: sum {u in units} buy[u] * reach[u];")
ampl$eval("subject to budget: sum {u in units} buy[u] * cost[u] <= p_budget;")
ampl$eval("subject to capacity: sum {u in units} buy[u] * person_weeks[u] <= max_person_weeks;")
ampl
$solve() ampl
HiGHS 1.7.0: HiGHS 1.7.0: optimal solution; objective 92
2 simplex iterations
0 barrier iterations
$getVariable("buy")$getValues() ampl
index0 buy.val
1 magazine_pages 20
2 tv_mins 40
$close() ampl
c. Adding channel radio
Model keeps the same as previous, but data changes.
# Adding channel to data
<- c("tv_mins", "magazine_pages", "radio_min")
units <- c(20000, 10000, 2000)
cost <- c(1.8, 1, 0.25)
reach <- c(10, 1, 1)
min_units <- c(1, 3, 1/7)
person_weeks
# loading class and setting solver
<- new(AMPL, env)
ampl $setOption("solver","HiGHS")
ampl
$eval("param p_budget := 1000000;")
ampl$eval("param max_person_weeks := 100;")
ampl# loading model
$read("models/1.1 advertising_campaigns.mod") # read model located in folder models
ampl
# setting data
$setData(
ampldata.frame(
units = units,
cost = cost,
reach = reach,
min_units = min_units,
person_weeks = person_weeks
1, "units"
),
) # formulation
$eval("var buy{u in units} >= min_units[u];")
ampl$eval("maximize audience: sum {u in units} buy[u] * reach[u];")
ampl$eval("subject to budget: sum {u in units} buy[u] * cost[u] <= p_budget;")
ampl$eval("subject to capacity: sum {u in units} buy[u] * person_weeks[u] <= max_person_weeks;")
ampl
$solve() ampl
HiGHS 1.7.0: HiGHS 1.7.0: optimal solution; objective 117.75
1 simplex iterations
0 barrier iterations
$getVariable("buy")$getValues() ampl
index0 buy.val
1 magazine_pages 1
2 radio_min 395
3 tv_mins 10
$close() ampl
2. Steel model with changes (amplpy)
Before solving execises, let’s review steel4 model explained in the chapter.
Steel4 model
from amplpy import AMPL, modules
# load all modules
modules.load() = AMPL() # instantiate AMPL object}
ampl "solver"] = "highs"
ampl.option["models/steel_original.mod") # read model located in folder models
ampl.read("data/steel_original.dat") # read dat located in folder models
ampl.read_data( ampl.solve()
HiGHS 1.7.0: HiGHS 1.7.0: optimal solution; objective 190071.4286
2 simplex iterations
0 barrier iterations
= ampl.getVariable("Make").getValues()
df print(df)
index0 | Make.val
'bands' | 3357.142857142858
'coils' | 500
'plate' | 3142.857142857143
ampl.reset() ampl.close()
Model developed in AMPL
```{AMPL}
set PROD; # products
set STAGE; # stages
param rate {PROD,STAGE} > 0; # tons per hour in each stage
param avail {STAGE} >= 0; # hours available/week in each stage
param profit {PROD}; # profit per ton
param commit {PROD} >= 0; # lower limit on tons sold in week
param market {PROD} >= 0; # upper limit on tons sold in week
var Make {p in PROD} >= commit[p], <= market[p]; # tons produced
maximize Total_Profit: sum {p in PROD} profit[p] * Make[p];
# Objective: total profits from all products
subject to Time {s in STAGE}:
sum {p in PROD} (1/rate[p,s]) * Make[p] <= avail[s];
# In each stage: total of hours used by all
# products may not exceed hours available
```
Data input in AMPL
```{AMPL}
set PROD := bands coils plate;
set STAGE := reheat roll;
param rate: reheat roll :=
bands 200 200
coils 200 140
plate 200 160 ;
param: profit commit market :=
bands 25 1000 6000
coils 30 500 4000
plate 29 750 3500 ;
param avail := reheat 35 roll 40 ;
```
a. Change constrain to equal
As the objective functions is a maximization, the problem would try to achive the equality as their is not another constraint. Therefore, changing the equality from <=
to =
does not change the result.
```{AMPL}
subject to Time {s in STAGE}:
sum {p in PROD} (1/rate[p,s]) * Make[p] = avail[s];
```
b. Max total weigth contraint
Adding a new constraint to put a upper bound to the variable Make.
# Loading class and setting solver
= AMPL() # instantiate AMPL object
ampl "solver"] = "highs"
ampl.option[
# Load model and data
"models/steel_original.mod") # read model located in folder models
ampl.read("data/steel_original.dat") # read dat located in folder models
ampl.read_data(
## Adding new constraint to limit tons manufactured.
eval("param max_weight := 6500;")
ampl.eval("subject to MaxWeight: sum {p in PROD} Make[p] <= max_weight;")
ampl.
# Solving problem
ampl.solve()
HiGHS 1.7.0: HiGHS 1.7.0: optimal solution; objective 183791.6667
3 simplex iterations
0 barrier iterations
= ampl.getVariable("Make").getValues()
df print(df)
index0 | Make.val
'bands' | 1541.666666666667
'coils' | 1458.333333333333
'plate' | 3500
ampl.close()
c. Produce as many tons as possible
The idea behind increase number of tons is to keep the objective functions, nevertheless, without the params profit.
```{AMPL}
maximize tons: sum {p in PROD} Make[p];
# Objective: total profits from all products
```
d. Lower bounds as constraints.
In this exercise we are asked to change the minimum of the tons make for each product, therefore, we required to redeclare the var to delete the lower bound.
# import libraries
import pandas as pd
# Loading class and setting solver
= AMPL() # instantiate AMPL object
ampl "solver"] = "highs"
ampl.option[
# Load model and data
"models/steel_original.mod") # read model located in folder models
ampl.read("data/steel_original.dat") # read dat located in folder models
ampl.read_data(
## Changing in formulation
eval("param shares {PROD};")
ampl.
# Creation of a list to input ampl
= {"bands": 0.4, "coils": 0.1, "plate": 0.4}
shares_list
# Inputing the list into the param shares created
"shares"] = {PROD: shares for PROD, (shares) in shares_list.items()}
ampl.param[
# Redeclare var Make to delete lower bound
eval("redeclare var Make {p in PROD} <= market[p];")
ampl.
# Create a new constraint with the lower bound by product
eval("subject to min_shares {j in PROD}: Make[j] >= shares[j] * sum {k in PROD} Make[k];")
ampl.
# Solving problem
ampl.solve()
HiGHS 1.7.0: HiGHS 1.7.0: optimal solution; objective 189700
5 simplex iterations
0 barrier iterations
= ampl.getVariable("Make").getValues()
df print(df)
index0 | Make.val
'bands' | 3500
'coils' | 700
'plate' | 2800
ampl.close()