Hello, World!

of mathematical programming

Introduction

Hello, World!

For any person who has ever taken a course or tutorial on a programming language, this last phrase is well known. Its tradition started with one of the first successful large-scale languages: C. In the book The C Programming Language (1978), authors Brian Kernighan and Dennis Ritchie selected the phrase to exemplify the program printing syntax, which the user of the program would be able to see at a terminal or screen display. With the evolution of programming languages, the need for simple examples for the software syntax persisted, so the Hello, World! became the go-to case of introductory programming.

Although mathematical modeling does not have a universal problem to present its concepts to beginners, there is a classical approach often used. It involves a small representation of a business decision, normally with a couple of constraints and levers, which can be easily solved by drawing them on a chart.

Inspired by those examples this ebook’s objective is to serve as an introduction to their combination in the use of mathematical programming languages. It shows how to translate a model of a mathematical problem for the following tools: Python, Java, C++, AMPL, and Excel.

The mathematical model

A manufacturer produces two types of wooden toys: soldiers and trains. A soldier sells for $27 and uses $10 worth of raw materials. Each soldier that is manufactured increases the variable labor and overhead costs by $14. A train sells for $21 and uses $9 worth of raw materials. Each train built increases the variable labor and overhead costs by $10. The manufacture of wooden soldiers and trains requires two types of skilled labor: carpentry and finishing. A soldier requires 2 hours of finishing labor and 1 hour of carpentry labor. A train requires 1 hour of finishing labor and 1 hour of carpentry labor.

Each week, the company can obtain all the needed raw materials, but only 100 finishing hours and 80 carpentry hours. Demand for trains is unlimited, but only 40 soldiers are bought, at the most, each week. The company wants to maximize weekly profit.

This description can be translated into the system of equations below.

Maximize Profit = (27-10-14)×Soldiers+ (21-9-10)×Trains

Subject to

2×Soldiers+Trains ≤100 (Finishing constraint)

Soldiers+Trains ≤80 (Carpentry constrain)

Soldiers ≤40 (Soldier's demand)

Soldiers ≥0 (Minimum number of soldiers)

Trains ≥0 (Minimum number of trains)

Implementation tools

Excel

Python

JAVA

C++

AMPL

Excel implementation

Although Excel is not a proper programming language, it is one of the most well known tools in the world and it has optimization options that can be used to solve our base case. The first thing to do is to make sure your Excel has the Solver Add-in selected to use. Then, the implementation of the problem begins by constructing tables for the inputs of the problem, as shown in Figure 2.

The next step consists of defining the outputs of the model, which are, in this case, dependent on the number of soldiers and trains built (Figure 3). Note that the file is using the naming cell feature, to ease the understanding of the formulas.

With both the input and output defined, it is possible to create the optimization model using the Excel wizard, which pops-up (Figure 4) after the Solver tool (Tools > Solver …) is activated. The first thing is to define the Set Objective field, which would be the H4, named Profit. Then, the To selection must be set as Max, since the model searches for the maximum possible profit. The third field By Changing Variable Cells wants to know the variables of the model, which are cells H9 and H10, named range Toy_production. After that, the Subject to the Constraints field needs to know each one of the problems bounds. Each constraint is added individually, like H7 (named Finishing_operation) (is lower or equal than) C15 (named Finishing), and so one. The last is to define the Select a Solving Method field, which for this linear problem is Simplex LP.

With the model fully completed, by calling for Solve, an optimal result is found (Figure 5).

Python implementation

Besides its great performance in data treatment, Python also has many libraries dedicated to optimization algorithms. This implementation uses Google's OR-Tools package, by importing

from ortools.linear_solver import pywraplp

After that, it is necessary to define the model function that will have all the model declarations.

def helloWorld():

The first command of the function is to define the object that encapsulates the model, including the solver (GLOP).

# Create the linear solver with the GLOP backend

solver = pywraplp.Solver.CreateSolver('GLOP')

The next step involves the declaration of the variables, using the function NumVar. To define each variable this function uses three arguments, in this order: minimum possible value, maximum possible value, and label. Therefore, the declaration is implicit in the definition of constraints.

# Create the variables soldiers and trains

soldiers = solver.NumVar(0, 40, 'Soldiers')

trains = solver.NumVar(0, 10000, 'Trains')

Followed by the declaration of constraints, use the functions Constraint and SetCoefficient. The former uses three arguments, in this order: minimum possible value, maximum possible value, and label; while the latter uses two arguments: the variable and its coefficient in the constraint.

# Create finishing constraint, 0 <= 2*soldiers + trains <= 100

finishing = solver.Constraint(0, 100, 'finishing')

finishing.SetCoefficient(soldiers, 2)

finishing.SetCoefficient(trains, 1)

# Create carpentry constraint, 0 <= soldiers + trains <= 80

carpentry = solver.Constraint(0, 80, 'carpentry')

carpentry.SetCoefficient(soldiers, 1)

carpentry.SetCoefficient(trains, 1)

Only the objective declaration is missing to complete the model, which is similar in format to the constraints

# Create the objective function, 3 * soldiers + 2 * trains

profit = solver.Objective()

profit.SetCoefficient(soldiers, 3)

profit.SetCoefficient(trains, 2)

profit.SetMaximization()

To run the model, call the Solve() function of the solver object.

# Solve the model

solver.Solve()

Finally, the results can be displayed in the interface by using the print command, like:

print('Profit = ', profit.Value())

print('Soldiers' = ', soldiers.solution_value.Value())

print('Trains' = ', trains.solution_value.Value())

Which returns the following

Profit = 180

Soldiers = 20

Trains = 60

JAVA implementation

Java is well known for being an object-oriented language that can run easily in all types of machines and strikes good balance between speed and learning curve. This implementation also uses Google's OR-Tools package by importing the following classes:

import com.google.ortools.Loader;

import com.google.ortools.linearsolver.MPConstraint;

import com.google.ortools.linearsolver.MPObjective;

import com.google.ortools.linearsolver.MPSolver;

import com.google.ortools.linearsolver.MPVariable;

Native libraries need loading with the following command

// Load libraries

Loader.loadNativeLibraries();

In typical java fashion, we will start by creating a Model class with a field for the solver, variables and objective function. We will include a constructor and a solve method in the class:

class Model {

private MPSolver solver; // LP solver

private MPVariable soldiers; // Soldiers variable

private MPVariable trains; // Trains variable

private MPObjective profit; // Objective function

public Model() {..}

public void solve() {...}

}

In the constructor, we first define the object that encapsulates the model, including the solver (GLOP).

// Initialize linear solver with GLOP backend

solver = MPSolver.createSolver("GLOP");

We also instantiate the variables and declare the constraints, which follow similar structures to the ones described in the Python section

// Create the variables soldiers and trains

soldiers = solver.makeNumVar(0, 40, "soldiers");

trains = solver.makeNumVar(0, 1000, "trains");

// Create finishing constraint, 0 <= 2*soldier + trains <= 100

MPConstraint finishing = solver.makeConstraint(0, 100, "finishing");

finishing.setCoefficient(soldiers, 2);

finishing.setCoefficient(trains, 1);

// Create carpentry constraint, 0 <= soldiers + trains <= 80

MPConstraint carpentry = solver.makeConstraint(0, 80, "carpentry");

carpentry.setCoefficient(soldiers, 1);

carpentry.setCoefficient(trains, 1);

To complete the model we instantiate the objective, which follows similar structures to the ones described in the Python section

// Create the objective function, 3 * soldiers + 2* trains

profit = solver.objective();

profit.setCoefficient(soldiers, 3);

profit.setCoefficient(trains, 2);

profit.setMaximization();

Finally we create a solve method that calls the OR tools solve and prints the result

public void solve() {

// Solve the model

solver.solve();

// Print objective function and solution

System.out.println("Solutions: ");

System.out.println("Objective value = " + profit.value());

System.out.println("soldiers = " + soldiers.solutionValue());

System.out.println("trains = " + trains.solutionValue());

}

Which returns the following

Solutions:

Objective value = 180

soldiers = 20

trains = 60

C++ implementation

The C and C++ languages are probably the ones used in the majority of the solver because they are high-performance languages. Similarly to the other implementations, we start by importing the significant libraries, and we also use the operations_research namespace

#include "ortools/linear_solver/linear_solver.h"

namespace operations_research {...}

Within the namespace created, we define our hello world function and declare the GLOP model

void helloWorld() {

// Create linear solver with GLOP backend

std::unique_ptr<MPSolver> solver(MPSolver::CreateSolver("GLOP"));

// …

}

We then declare our variables and constraints, which follow structures that are similar to the ones in the previous sections.

// Create the variables soldiers and trains

MPVariable* const soldiers = solver->MakeNumVar(0, 40, "Soldiers");

MPVariable* const trains = solver->MakeNumVar(0, 1000, "Trains");

// Create finishing constraint, 0 <= 2*soldier + trains <= 100

MPConstraint* const finishing = solver->MakeRowConstraint(0,100,"finishing");

finishing->SetCoefficient(soldiers, 2);

finishing->SetCoefficient(trains, 1);

// Create carpentry constraint, 0 <= soldiers + trains <= 80

MPConstraint* const carpentry = solver->MakeRowConstraint(0, 80, "carpentry");

carpentry->SetCoefficient(soldiers, 1);

carpentry->SetCoefficient(trains, 1);

To complete the model we declare the objective function, which follow structures that are similar to the ones in the previous sections.

// Create the objective function, 3 * soldiers + 2* trains

MPObjective* const objective = solver->MutableObjective();

objective->SetCoefficient(soldiers, 3);

objective->SetCoefficient(trains, 2);

objective->SetMaxiization();

To run the model, call the Solve() function of the solver object.

// Solve the model

solver->Solve();

std::cout << "Solutions:" << std::endl;

std::cout << "Objective value = " << objective->Value() << std::endl;

std::cout << "soldiers = " << soldiers->solution_value() << std::endl;

std::cout << "trains = " << trains->solution_value() << std::endl;

Finally, we create a main function of the program, and call for the helloWorld() function.

int main() {

operations_research::helloWorld();

return 0;

}

Which returns the following

Solutions:

Objective value = 180

soldiers = 20

trains = 60

AMPL implementation

The A Mathematical Programming Language (AMPL) is dedicated to algebraic representation of models and it is very similar to the formulation of the problem. Once this software is installed in the machine, open its Terminal and open AMPL using the command ampl. From there, the declaration of the model is straightforward. The first thing to do is to define the sets of the model:

set Toys;

set Costs;

set Operations;

Then the parameters:

param production_costs {t in Toys, c in Costs};

param price {t in Toys};

param demand {t in Toys};

param operating_time {t in Toys, o in Operations};

param working_hours {o in Operations};

Followed by variables:

var Production {t in Toys} >= 0;

Then the model's constraints:

subject to Operations_constraints {o in Operations}:

sum{t in Toys} (Production[t] * operating_time[t,o]) <= working_hours[o];

subject to Demand_constraints {t in Toys}: Production[t] <= demand[t];

And, finally, the objective function:

maximize Profit: sum{t in Toys} (Production[t] * (price[t] - sum{c in Costs} production_costs[t,c]));

With all structures in place, the only thing left is to define the elements of sets and value of parameters, using the command data:

data;

set Toys := 'Soldier' 'Trains';

set Costs := 'Raw materials' 'Labor and overhead';

set Operations := 'Finishing' 'Carpentry';

param production_costs :

'Raw materials' 'Labor and overhead' :=

'Soldier' 10 14

'Trains' 9 10;

param: price demand :=

'Soldier' 27 40

'Trains' 21 10000;

param operating_time :

'Finishing' 'Carpentry' :=

'Soldier' 2 1

'Trains' 1 1;

param working_hours :=

'Finishing' 100

'Carpentry' 80;

Finally, the programmer selects a solver, using the command option, and finds a solution, by calling the command solve. For instance, we use Xpress as the solver engine.

option solver xpress;

solve;

Which prints the solution in the Terminal:

XPRESS 8.5(33.01.05): Optimal solution found

Objective 180

2 simplex iterations

Technology comparison

Although all approaches delivered the same results for the simple case, it may not be the case for other models. The technological choice for development must take into account several more factors, such as: model complexity, processing capacity, programming skills, reserved budget, graphical user interface (GUI), versioning system, etc.

For instance, Excel is a great option for simple optimization problems, being able to easily generate visualizations of the solution and change parameters’ values. However, it has a low number of variables it can consider and low processing capacity, so larger problems are usually out of scope. Additionally, it does not have any versioning system, which can be confusing if many coworkers should use the same model.

The three main programming languages listed, Python, Java, and C/C++, are very similar in their capabilities to solve problems, and are able to solve problems of high complexity and with powerful processing power. In their advantage, they can be used for free, and are compatible with versioning technologies, like git. However, creating an application to show results and interact with users demands excellent programming skills.

The main benefit of the AMPL software is its clear division of model and data, but it is the only option that requires a license, which is free for students and paid for researchers and companies. AMPL enables non-programmers to easily model complex problems and still be able to reach the same performance levels and versioning systems of other programming languages. Additionally it can count on a flexible interface with charts and tables, called Quandec, which has an in-built data structure for scenario management, online access, parallel collaboration, comparison tools, and much more.

Conclusion

The world is filled with optimization opportunities that just need the right tool for the job. We hope this e-book has added another one into your toolbox. However, if the challenge seems to demand more experience in mathematical programming, please reach out to Cassotis at our channels and we can start a collaboration.