Note: This page has been translated by MathWorks. Click here to see

To view all translated materials including this page, select Country from the country navigator on the bottom of this page.

To view all translated materials including this page, select Country from the country navigator on the bottom of this page.

This example shows how to convert a linear problem from mathematical form into Optimization Toolbox™ solver syntax using the problem-based approach.

The variables and expressions in the problem represent a model of operating a chemical plant, from an example in Edgar and Himmelblau [1]. There are two videos that describe the problem.

Mathematical Modeling with Optimization, Part 1 shows the problem in pictorial form. It shows how to generate the mathematical expressions of Model Description from the picture.

Optimization Modeling, Part 2: Problem-Based Solution of a Mathematical Model describes how to convert these mathematical expressions into Optimization Toolbox solver syntax. This video shows how to solve the problem, and how to interpret the results.

The remainder of this example is concerned solely with transforming the problem to solver syntax. The example closely follows the video Optimization Modeling, Part 2: Problem-Based Solution of a Mathematical Model.

The video Mathematical Modeling with Optimization, Part 1 suggests that one way to convert a problem into mathematical form is to:

Get an overall idea of the problem

Identify the goal (maximizing or minimizing something)

Identify (name) variables

Identify constraints

Determine which variables you can control

Specify all quantities in mathematical notation

Check the model for completeness and correctness

For the meaning of the variables in this section, see the video Mathematical Modeling with Optimization, Part 1.

The optimization problem is to minimize the objective function, subject to all the other expressions as constraints.

The objective function is:

`0.002614 HPS + 0.0239 PP + 0.009825 EP`

.

The constraints are:

`2500`

≤ `P1`

≤
`6250`

`I1`

≤
`192,000`

`C`

≤
`62,000`

`I1 - HE1`

≤
`132,000`

```
I1 = LE1 + HE1 +
C
```

```
1359.8 I1 = 1267.8 HE1 + 1251.4 LE1 +
192 C + 3413 P1
```

`3000`

≤
`P2`

≤
`9000`

`I2`

≤
`244,000`

`LE2`

≤
`142,000`

```
I2 = LE2 +
HE2
```

```
1359.8 I2 = 1267.8 HE2 + 1251.4 LE2
+ 3413 P2
```

```
HPS = I1 + I2 +
BF1
```

```
HPS = C + MPS +
LPS
```

```
LPS = LE1 + LE2 +
BF2
```

```
MPS = HE1 + HE2 + BF1 -
BF2
```

`P1 + P2 + PP`

≥
`24,550`

`EP + PP`

≥
`12,000`

`MPS`

≥
`271,536`

`LPS`

≥
`100,623`

All variables are
positive.

The first solution method involves creating an optimization variable for each problem variable. As you create the variables, include their bounds.

P1 = optimvar('P1','LowerBound',2500,'UpperBound',6250); P2 = optimvar('P2','LowerBound',3000,'UpperBound',9000); I1 = optimvar('I1','LowerBound',0,'UpperBound',192000); I2 = optimvar('I2','LowerBound',0,'UpperBound',244000); C = optimvar('C','LowerBound',0,'UpperBound',62000); LE1 = optimvar('LE1','LowerBound',0); LE2 = optimvar('LE2','LowerBound',0,'UpperBound',142000); HE1 = optimvar('HE1','LowerBound',0); HE2 = optimvar('HE2','LowerBound',0); HPS = optimvar('HPS','LowerBound',0); MPS = optimvar('MPS','LowerBound',271536); LPS = optimvar('LPS','LowerBound',100623); BF1 = optimvar('BF1','LowerBound',0); BF2 = optimvar('BF2','LowerBound',0); EP = optimvar('EP','LowerBound',0); PP = optimvar('PP','LowerBound',0);

Create an optimization problem container. Include the objective function in the problem.

`linprob = optimproblem('Objective',0.002614*HPS + 0.0239*PP + 0.009825*EP);`

There are three linear inequalities in the problem expressions:

`I1 - HE1`

≤
`132,000`

`EP + PP`

≥
`12,000`

`P1 + P2 + PP`

≥
`24,550`

.

Create these inequality constraints and include them in the problem.

linprob.Constraints.cons1 = I1 - HE1 <= 132000; linprob.Constraints.cons2 = EP + PP >= 12000; linprob.Constraints.cons3 = P1 + P2 + PP >= 24550;

There are eight linear equalities:

`I2 = LE2 + HE2`

```
LPS = LE1 +
LE2 + BF2
```

```
HPS = I1 + I2 +
BF1
```

```
HPS = C + MPS +
LPS
```

```
I1 = LE1 + HE1 +
C
```

```
MPS = HE1 + HE2 + BF1 -
BF2
```

```
1359.8 I1 = 1267.8 HE1 + 1251.4 LE1
+ 192 C + 3413 P1
```

```
1359.8 I2 = 1267.8 HE2
+ 1251.4 LE2 + 3413 P2
```

.

Include these constraints as well.

linprob.Constraints.econs1 = LE2 + HE2 == I2; linprob.Constraints.econs2 = LE1 + LE2 + BF2 == LPS; linprob.Constraints.econs3 = I1 + I2 + BF1 == HPS; linprob.Constraints.econs4 = C + MPS + LPS == HPS; linprob.Constraints.econs5 = LE1 + HE1 + C == I1; linprob.Constraints.econs6 = HE1 + HE2 + BF1 == BF2 + MPS; linprob.Constraints.econs7 = 1267.8*HE1 + 1251.4*LE1 + 192*C + 3413*P1 == 1359.8*I1; linprob.Constraints.econs8 = 1267.8*HE2 + 1251.4*LE2 + 3413*P2 == 1359.8*I2;

The problem formulation is complete. Solve the problem using
`solve`

.

linsol = solve(linprob);

Optimal solution found.

Evaluate the objective function. (You could have asked for this value when you
called `solve`

.)

evaluate(linprob.Objective,linsol)

ans = 1.2703e+03

The lowest-cost method of operating the plant costs $1,207.30.

Examine the solution variable values.

tbl = struct2table(linsol)

tbl = 1×16 table BF1 BF2 C EP HE1 HE2 HPS I1 I2 LE1 LE2 LPS MPS P1 P2 PP ___ ___ ______ ______ __________ __________ __________ __________ ________ ___ __________ __________ __________ ____ ______ _____ 0 0 8169.7 760.71 1.2816e+05 1.4338e+05 3.8033e+05 1.3633e+05 2.44e+05 0 1.0062e+05 1.0062e+05 2.7154e+05 6250 7060.7 11239

This table is too wide to see easily. Stack the variables to get them to a vertical orientation.

vars = {'P1','P2','I1','I2','C','LE1','LE2','HE1','HE2',... 'HPS','MPS','LPS','BF1','BF2','EP','PP'}; outputvars = stack(tbl,vars,'NewDataVariableName','Amt','IndexVariableName','Var')

outputvars = 16×2 table Var Amt ___ __________ P1 6250 P2 7060.7 I1 1.3633e+05 I2 2.44e+05 C 8169.7 LE1 0 LE2 1.0062e+05 HE1 1.2816e+05 HE2 1.4338e+05 HPS 3.8033e+05 MPS 2.7154e+05 LPS 1.0062e+05 BF1 0 BF2 0 EP 760.71 PP 11239

`BF1`

,`BF2`

, and`LE1`

are`0`

, their lower bounds.`I2`

is`244,000`

, its upper bound.The nonzero components of the objective function (cost) are

`HPS`

—`380,328.74`

`PP`

—`11,239.29`

`EP`

—`760.71`

The video Optimization Modeling, Part 2: Problem-Based Solution of a Mathematical Model gives interpretations of these characteristics in terms of the original problem.

Alternatively, you can solve the problem using just one optimization variable that has indices with the names of the problem variables. This method enables you to give a lower bound of zero to all problem variables at once.

vars = {'P1','P2','I1','I2','C','LE1','LE2','HE1','HE2',... 'HPS','MPS','LPS','BF1','BF2','EP','PP'}; x = optimvar('x',vars,'LowerBound',0);

Include the bounds on the variables using dot notation.

x('P1').LowerBound = 2500; x('P2').LowerBound = 3000; x('MPS').LowerBound = 271536; x('LPS').LowerBound = 100623; x('P1').UpperBound = 6250; x('P2').UpperBound = 9000; x('I1').UpperBound = 192000; x('I2').UpperBound = 244000; x('C').UpperBound = 62000; x('LE2').UpperBound = 142000;

The remainder of the problem setup is similar to the setup using separate
variables. The difference is that, instead of addressing a variable by its name,
such as `P1`

, you address it using its index,
`x('P1')`

.

Create the problem object, include the linear constraints, and solve the problem.

linprob = optimproblem('Objective',0.002614*x('HPS') + 0.0239*x('PP') + 0.009825*x('EP')); linprob.Constraints.cons1 = x('I1') - x('HE1') <= 132000; linprob.Constraints.cons2 = x('EP') + x('PP') >= 12000; linprob.Constraints.cons3 = x('P1') + x('P2') + x('PP') >= 24550; linprob.Constraints.econs1 = x('LE2') + x('HE2') == x('I2'); linprob.Constraints.econs2 = x('LE1') + x('LE2') + x('BF2') == x('LPS'); linprob.Constraints.econs3 = x('I1') + x('I2') + x('BF1') == x('HPS'); linprob.Constraints.econs4 = x('C') + x('MPS') + x('LPS') == x('HPS'); linprob.Constraints.econs5 = x('LE1') + x('HE1') + x('C') == x('I1'); linprob.Constraints.econs6 = x('HE1') + x('HE2') + x('BF1') == x('BF2') + x('MPS'); linprob.Constraints.econs7 = 1267.8*x('HE1') + 1251.4*x('LE1') + 192*x('C') + 3413*x('P1') == 1359.8*x('I1'); linprob.Constraints.econs8 = 1267.8*x('HE2') + 1251.4*x('LE2') + 3413*x('P2') == 1359.8*x('I2'); [linsol,fval] = solve(linprob);

Optimal solution found.

Examine the solution as a vertical table.

tbl = table(vars',linsol.x')

tbl = 16×2 table Var1 Var2 _____ __________ 'P1' 6250 'P2' 7060.7 'I1' 1.3633e+05 'I2' 2.44e+05 'C' 8169.7 'LE1' 0 'LE2' 1.0062e+05 'HE1' 1.2816e+05 'HE2' 1.4338e+05 'HPS' 3.8033e+05 'MPS' 2.7154e+05 'LPS' 1.0062e+05 'BF1' 0 'BF2' 0 'EP' 760.71 'PP' 11239

[1] Edgar, Thomas F., and
David M. Himmelblau. *Optimization of Chemical Processes.*
McGraw-Hill, New York, 1988.