({{, (RefLink[LinearOptimization, paclet:ref/LinearOptimization])[f,cons,vars] 
finds values of variables vars that minimize the linear objective f subject to linear constraints cons.}, {, (RefLink[LinearOptimization, paclet:ref/LinearOptimization])[c,{a,b}] 
finds a vector x that minimizes the linear objective c.x subject to the linear inequality constraints a.x+b⪰0.}, {, (RefLink[LinearOptimization, paclet:ref/LinearOptimization])[c,{a,b},{(Subscript[a, eq]),(Subscript[b, eq])}]
includes the linear equality constraints (Subscript[a, eq]).x+(Subscript[b, eq])==0.}, {, (RefLink[LinearOptimization, paclet:ref/LinearOptimization])[…,"prop"]
specifies what solution property "prop" should be returned.}})

OpenersNotesSection(Spacer1)Details and Options

Linear optimization is also known as linear programming (LP).

Linear optimization is a convex optimization problem that can be solved globally and efficiently.

Linear optimization finds x that solves the primal problem: »

The constraints cons can be specified by:

With (RefLink[LinearOptimization, paclet:ref/LinearOptimization])[f,cons,vars], parameter equations of the form par==val, where par is not in vars and val is numerical or an array with numerical values, may be included in the constraints to define parameters used in f or cons.

The primal minimization problem has a related maximization problem that is the Lagrangian dual problem. The dual maximum value is always less than or equal to the primal minimum value, so it provides a lower bound. The dual maximizer provides information about the primal problem, including sensitivity of the minimum value to changes in the constraints.

The Lagrangian dual problem for linear optimization is given by: »

For linear optimization, strong duality always holds, meaning that if there is a solution to the primal minimization problem, then there is a solution to the dual maximization problem, and the dual maximum value is equal to the primal minimum value.

The possible solution properties "prop" include:

The following options can be given:

The option (RefLink[Method, paclet:ref/Method])->method may be used to specify the method to use. Available methods include:

With (RefLink[WorkingPrecision, paclet:ref/WorkingPrecision])->(RefLink[Automatic, paclet:ref/Automatic]), the precision is taken automatically from the precision of the input arguments unless a method is specified that only works with machine precision, in which case machine precision is used.

OpenersPrimaryExamplesSection(Spacer1)Examples  (74)

OpenersExampleSection0(Spacer1)Basic Examples  (2)

Minimize x+y subject to the constraints x+2y≥3, x≥0, y≥0:

res=LinearOptimization[x+y,{x+2y≥3,x≥0,y≥0},{x,y}]

The optimal point lies in a region defined by the constraints and where x+y is smallest within the region:

Show[Plot3D[x+y,{x,y}∈ImplicitRegion[x+2 y≥3&&x≥0&&y≥0,{x,y}],],Graphics3D[{Blue, PointSize[0.04],Point[{x,y,x+y}/.res]}]]

Minimize x+y subject to the constraints x+2y==3,x≥0,y≥0:

res=LinearOptimization[x+y,{x+2y==3,x≥0,y≥0},{x,y}]

Use the equivalent matrix representation:

LinearOptimization[{1,1},{{{1,0},{0,1}},{0,0}},{{{1,2}},{-3}}]

OpenersExampleSection1(Spacer1)Scope  (26)

OpenersExampleSubsection0(Spacer1)Basic Uses  (16)

Minimize x+y subject to the constraints x+2y≥3, x≥-1,y≥-1:

res=LinearOptimization[x+y,{x+2y≥3,x≥-1,y≥-1},{x,y}]

Minimize x-y subject to the constraints x+y+z=1/2,x-2z=1 and 2x-y≥1:

LinearOptimization[x-y,{x+y+z==1/2,x-2z==1,2x-y≥1},{x,y,z}]

Minimize x+y subject to the constraints x+2y≥3 and the bounds -1≤x≤2, -1≤y≤1:

LinearOptimization[x+y, {x+2y≥3,-1≤x≤2,-1≤y≤1},{x,y}]

Use (RefLink[VectorLessEqual, paclet:ref/VectorLessEqual]) to express several (RefLink[LessEqual, paclet:ref/LessEqual]) inequality constraints at once:

LinearOptimization[x+y,VectorLessEqual[{{x-2y,-x+y},{3,2}}],{x,y}]

Use (Key0[esc]) v<= (Key0[esc]) to enter the vector inequality in a compact form:

LinearOptimization[x+y,{x-2y,-x+y}\[VectorLessEqual]{3,2},{x,y}]

An equivalent form using scalar inequalities:

LinearOptimization[x+y,{x-2y≤3,-x+y≤2},{x,y}]

Use (RefLink[VectorGreaterEqual, paclet:ref/VectorGreaterEqual]) to express several (RefLink[GreaterEqual, paclet:ref/GreaterEqual]) inequality constraints at once:

LinearOptimization[x+y,VectorGreaterEqual[{{x+2y,x},{3,-1}}],{x,y}]

Use (Key0[esc]) v>= (Key0[esc]) to enter the vector inequality in a compact form:

LinearOptimization[x+y,{x+2y,x}\[VectorGreaterEqual]{3,-1},{x,y}]

Use (RefLink[Equal, paclet:ref/Equal]) to express several equality constraints at once:

LinearOptimization[x-y,{x+y,x-y}=={1/2,1},{x,y}]

An equivalent form using several scalar equalities:

LinearOptimization[x-y,{x+y==1/2,x-y==1},{x,y}]

Specify the constraints using a combination of scalar and vector inequalities:

LinearOptimization[x+y,{x+2y≥3,-1\[VectorLessEqual]{x,y}\[VectorLessEqual]{2,1}},{x,y}]

An equivalent form using scalar inequalities:

LinearOptimization[x+y,{x+2y≥3,-1≤x≤2,-1≤y≤1},{x,y}]

Minimize 3x+4y subject to x-2y≤3,-x+y≤2. Use a vector variable v={x,y}:

LinearOptimization[{3,4}.v,{{1,-2}.v≤3,{-1,1}.v≤2},v]

An equivalent form using scalar variables:

LinearOptimization[{3,4}.{x,y},{{1,-2}.{x,y}≤3,{-1,1}.{x,y}≤2},{x,y}]

Use vector variable and vector inequalities to specify the problem:

LinearOptimization[{3,4}.v,{{1,-2},{-1,1}}.v\[VectorLessEqual]{3,2},v]

Minimize the function (∑_(i=1)^3)(Subscript[x, i]) subject to the constraint a.x\[VectorGreaterEqual]b, (Subscript[x, 1])≥1,x∈(ℝ^3):

{a,b}={{{1,1,0},{0,1,1},{1,1,1}},{2,2,-1}};

Use (RefLink[Indexed, paclet:ref/Indexed]) to access components of a vector variable, e.g. (IndexedDefault[x, 1]):

LinearOptimization[Total[x],{a.x\[VectorGreaterEqual]b,Indexed[x,1]≥1},x]

Use parameter equations to specify the coefficients of the objective and constraints:

parameterEqs={c=={3,4},a1=={1,-2},a2=={-1,1}};
constraints={a1.{x,y}≤3,a2.{x,y}≤2};

LinearOptimization[c.{x,y},{parameterEqs,constraints},{x,y}]

Use parameter equations to specify coefficients of several constraints:

parameterEqs={c=={3,4},a=={{1,-2},{-1,1}},b=={3,2}};

LinearOptimization[c.{x,y},{parameterEqs,a.{x,y}\[VectorLessEqual]b},{x,y}]

Using vector variables:

LinearOptimization[c.v,{parameterEqs,a.v\[VectorLessEqual]b},v]

Minimize the function (∑_(i=1)^2)(Subscript[x, i])+(∑_(j=1)^3)(Subscript[y, i]) subject to a.x\[VectorGreaterEqual]b, y\[VectorGreaterEqual]1, y∈(ℝ^3):

parameterEqs={a=={{1,2},{1,0}},b=={3,-1}};

LinearOptimization[Total[x]+Total[y],{parameterEqs,a.x\[VectorGreaterEqual]b,y\[VectorGreaterEqual]1},{x,y}]

Use (RefLink[Vectors, paclet:ref/Vectors])[n] to specify the dimension of vector variables when needed:

LinearOptimization[Total[x]+Total[y],{parameterEqs,a.x\[VectorGreaterEqual]b,y\[VectorGreaterEqual]1},{x,y∈Vectors[3]}]

Specify non-negative constraints using (RefLink[NonNegativeReals, paclet:ref/NonNegativeReals]) ((NonNegativeReals)):

LinearOptimization[Total[v],{1,2}.v≥3,v∈Vectors[2,NonNegativeReals]]

An equivalent form using vector inequalities:

LinearOptimization[Total[v],{{1,2}.v≥3,v\[VectorGreaterEqual]0},v]

Specify non-positive constraints using (RefLink[NonPositiveReals, paclet:ref/NonPositiveReals]) ((NonPositiveReals)):

LinearOptimization[-Total[v],{1,2}.v\[VectorLessEqual]-3,v∈Vectors[2,NonPositiveReals]]

An equivalent form using vector inequalities:

LinearOptimization[-Total[v],{{1,2}.v\[VectorLessEqual]-3,v\[VectorLessEqual]0},v]

Minimize c.x subject to a.x+b\[VectorGreaterEqual]0 by directly giving the matrix a and vectors b and c:

{a,b,c}={{{1,2},{1,0},{0,1}}, {-3,1,1},{1,1}};

LinearOptimization[c,{a,b}]

The equivalent form using variables and inequalities:

constraints=Thread[a.{x,y}+b≥{0,0,0}]

LinearOptimization[c.{x,y},constraints,{x,y}]

OpenersExampleSubsection1(Spacer1)Primal Model Properties  (3)

Minimize the function x+y subject to the constraints x+2y≥3, -1≤x≤2,-1≤y≤1:

LinearOptimization[x+y,{x+2y>=3,-1≤x≤2,-1≤y≤1},{x,y}]

Get the primal minimizer as a vector:

LinearOptimization[x+y,{x+2y>=3,-1≤x≤2,-1≤y≤1},{x,y}, "PrimalMinimizer"]

Get the minimal value:

LinearOptimization[x+y,{x+2y>=3,-1≤x≤2,-1≤y≤1},{x,y}, "PrimalMinimumValue"]

Obtain the primal minimum value using symbolic inputs:

obj=x+y+1;
cons= {x+2y≥3,-1≤x≤2,-1≤y≤1};

(Subscript[p, sym])=LinearOptimization[obj,cons,{x,y},"PrimalMinimumValue"]

Extract the matrix-vector inputs of the optimization problem:

{c,{a,b}}=LinearOptimization[obj,cons,{x,y},{"ObjectiveVector","LinearInequalityConstraints"}];

Solve using the matrix-vector form:

(Subscript[p, m])=LinearOptimization[c,{a,b},"PrimalMinimumValue"]

Recover the symbolic primal value by adding the objective constant:

const=LinearOptimization[obj,cons,{x,y},"ObjectiveConstant"]

(Subscript[p, sym])==(Subscript[p, m])+const

Find the slack associated with a minimization problem:

obj=x+y;
cons={x+2y==3,-1≤x≤2,-1≤y≤1};

{(Subscript[s, ineq]),(Subscript[s, eq])}=LinearOptimization[obj,cons,{x,y},"Slack"]

Get the primal minimizer and the constraint matrices:

{(SuperStar[x]),{a,b},{(Subscript[a, eq]),(Subscript[b, eq])}}=LinearOptimization[obj,cons,{x,y},];

The slack for inequality constraints a.(SuperStar[x])+b\[VectorGreaterEqual]0 is a vector (Subscript[s, ineq]) such that a.(SuperStar[x])+b-(Subscript[s, ineq])==0:

Norm[a.(SuperStar[x])+b-(Subscript[s, ineq])]==0

The slack for the equality constraints a.(SuperStar[x])+b==0 is a vector (Subscript[s, eq]) such that a.(SuperStar[x])+b-(Subscript[s, eq])==0:

Norm[(Subscript[a, eq]).(SuperStar[x])+(Subscript[b, eq])-(Subscript[s, eq])]==0

The equality slack (Subscript[s, eq]) is typically a zero vector:

(Subscript[s, eq])

OpenersExampleSubsection2(Spacer1)Dual Model Properties  (3)

Minimize c.x subject to a.x+b\[VectorGreaterEqual]0 and (Subscript[a, eq]).x+(Subscript[b, eq])==0:

c={5,2,1};
{a,b}={{{1,1,0},{1,-1,0}},{-1,1}};
{(Subscript[a, eq]),(Subscript[b, eq])}={{{1,2,0},{1,1,1}},{-3,2}};
xv=Table[Indexed[x,i],{i,3}];

psol=LinearOptimization[c.xv,{a.xv+b\[VectorGreaterEqual]0,(Subscript[a, eq]).xv+(Subscript[b, eq])==0},xv]

The dual problem is to maximize -(Subscript[b, eq]).ν-b.λ subject to (Transpose[((Subscript[a, eq]))]).ν+(Transpose[a]).λ==c and λ\[VectorGreaterEqual]0:

{λv,νv}={Table[Indexed[λ,i],{i,2}],Table[Indexed[ν,i],{i,2}]};

dsol=LinearOptimization[b.λv+(Subscript[b, eq]).νv,{a\[Transpose].λv+(Subscript[a, eq])\[Transpose].νv==c,λv\[VectorGreaterEqual]0},Join[λv,νv]]

The primal minimum value and the dual maximum value coincide because of strong duality:

{c.xv/.psol,-(b.λv+(Subscript[b, eq]).νv)/.dsol}

That is the same as having a duality gap of zero. In general, c.(SuperStar[x])≤-b.(SuperStar[λ])-(Subscript[b, eq]).(SuperStar[ν]) at optimal points:

LinearOptimization[c.xv,{a.xv+b\[VectorGreaterEqual]0,(Subscript[a, eq]).xv+(Subscript[b, eq])==0},xv,"DualityGap"]

Construct the dual problem using constraint matrices extracted from the primal problem:

obj=5x+2y+z;
cons={x+y≥1,x-y≥-1,x+2y==3,x+y+z==-2};

LinearOptimization[obj,cons,{x,y,z}]

Extract the objective vector and constraint matrices:

{c,{a,b},{(Subscript[a, eq]),(Subscript[b, eq])}}=LinearOptimization[obj,cons,{x,y,z},{"ObjectiveVector","LinearInequalityConstraints","LinearEqualityConstraints"}];

The dual problem is to maximize -(Subscript[b, eq]).ν-b.λ subject to (Transpose[((Subscript[a, eq]))]).ν+(Transpose[a]).λ==c and λ\[VectorGreaterEqual]0:

LinearOptimization[b.λ+(Subscript[b, eq]).ν,{a\[Transpose].λ+(Subscript[a, eq])\[Transpose].ν==c,λ\[VectorGreaterEqual]0},{λ,ν}]

Get the dual maximum value directly using solution properties:

LinearOptimization[5x+2y+z,{x+y≥1,x-y≥-1,x+2y==3,x+y+z==-2},{x,y,z},"DualMaximumValue"]

Get the dual maximizer directly using solution properties:

LinearOptimization[5x+2y+z,{x+y≥1,x-y≥-1,x+2y==3,x+y+z==-2},{x,y,z},"DualMaximizer"]

OpenersExampleSubsection3(Spacer1)Sensitivity Properties  (4)

Use "ConstraintSensitivity" to find the change in optimal value due to constraint relaxations:

c={1,-1};
{a,b}={{{-1,-2},{-2,-1},{1,0},{0,1}},{2.,2.,0,0}};
(SuperStar[p])=LinearOptimization[c,{a,b},"PrimalMinimumValue"];

The first vector is inequality sensitivity and the second is equality sensitivity:

{(Subscript[s, λ]),(Subscript[s, ν])}=LinearOptimization[c,{a,b},"ConstraintSensitivity"]

Consider new constraints a.x+b+δb\[VectorGreaterEqual]0 where δb is the relaxation:

δb={1,0,1/2,0};

The expected new optimal value is given by:

(SuperStar[p])+δb.(Subscript[s, λ])

Compare to directly solving the relaxed problem:

LinearOptimization[c,{a,b+δb},"PrimalMinimumValue"]

Each sensitivity is associated with an inequality or equality constraint:

{(Subscript[s, λ]),(Subscript[s, ν])}=LinearOptimization[x-y+z,{x+2y≤2,2x+y≤2,x≥0,y≥0,
x+y+z==1},{x,y,z},"ConstraintSensitivity"]

Extract the constraints:

{{a,b},{(Subscript[a, eq]),(Subscript[b, eq])}}=LinearOptimization[x-y+z,{x+2y≤2,2x+y≤2,x≥0,y≥0,x+y+z==1},{x,y,z},{"LinearInequalityConstraints","LinearEqualityConstraints"}];

The inequality constraints and their associated sensitivity:

Table[{a[[i]].{x,y,z}+b[[i]]≥0,(Subscript[s, λ])[[i]]},{i,4}]

The equality constraints and their associated sensitivity:

Table[{(Subscript[a, eq])[[i]].{x,y,z}+(Subscript[b, eq])[[i]]==0,(Subscript[s, ν])[[i]]},{i,1}]

The change in optimal value due to constraint relaxation is proportional to the sensitivity value:

c={1,-1,1};
{a,b}={{{2,-1,0},{1,1,0},{0,1,0}},{1,0,0}};
{(Subscript[a, eq]),(Subscript[b, eq])}={{{1,1,1},{1,0,-2}},-{1/2,1}};

Compute the minimal value and constraint sensitivity:

{(SuperStar[p]),{(Subscript[s, λ]),(Subscript[s, ν])}}=LinearOptimization[c,{a,b},{(Subscript[a, eq]),(Subscript[b, eq])}, {"PrimalMinimumValue","ConstraintSensitivity"}]

A zero sensitivity will not change the optimal value if the constraint is relaxed:

(Subscript[s, λ])[[2]]

(SuperStar[p])==LinearOptimization[c,{a,b+{0,1,0}},{(Subscript[a, eq]),(Subscript[b, eq])},"PrimalMinimumValue"]

A negative sensitivity will decrease the optimal value:

(Subscript[s, λ])[[1]]

(SuperStar[p])>LinearOptimization[c,{a,b+{1,0,0}},{(Subscript[a, eq]),(Subscript[b, eq])},"PrimalMinimumValue"]

A positive sensitivity will increase the optimal value:

(Subscript[s, ν])[[1]]

(SuperStar[p])<LinearOptimization[c,{a,b},{(Subscript[a, eq]),(Subscript[b, eq])+{1,0}},"PrimalMinimumValue"]

The "ConstraintSensitivity" is related to the dual maximizer of the problem:

{{(Subscript[s, λ]),(Subscript[s, ν])},{(SuperStar[λ]),(SuperStar[ν])}}=LinearOptimization[x-y+z,{x+2y≤2,2x+y≤2,x≥0,y≥0,
x+y+z==1},{x,y,z},{"ConstraintSensitivity","DualMaximizer"}]

The inequality sensitivity (Subscript[s, λ])==((∂(SuperStar[p]))/(∂b)) satisfies (Subscript[s, λ])==-(SuperStar[λ]), where (SuperStar[λ]) is the dual inequality maximizer:

(Subscript[s, λ])==-(SuperStar[λ])

The equality sensitivity (Subscript[s, ν])==((∂(SuperStar[p]))/(∂(Subscript[b, eq]))) satisfies (Subscript[s, ν])==-(SuperStar[ν]), where (SuperStar[ν]) is the dual equality maximizer:

(Subscript[s, ν])==-(SuperStar[ν])

OpenersExampleSection2(Spacer1)Options  (9)

OpenersExampleSubsection4(Spacer1)Method  (5)

The default method for (RefLink[MachinePrecision, paclet:ref/MachinePrecision]) input is "CLP":

obj=x+y;
constr={x+2 y≥3,-1≤x≤2,-1≤y≤1};

LinearOptimization[obj,N[constr],{x,y},Method->Automatic]

LinearOptimization[obj,N[constr],{x,y},Method->"CLP"]

The default method for arbitrary or infinite (RefLink[WorkingPrecision, paclet:ref/WorkingPrecision]) is "Simplex":

LinearOptimization[obj,constr,{x,y},Method->Automatic]

LinearOptimization[obj,constr,{x,y},Method->"Simplex"]

All methods work for (RefLink[MachinePrecision, paclet:ref/MachinePrecision]) input:

obj=x+y;
constr={x+2 y≥3,-1≤x≤2,-1≤y≤1};

Table[LinearOptimization[obj, N[constr],{x,y},Method->met],{met,{"Simplex","RevisedSimplex","InteriorPoint","CLP"}}]

"Simplex" and "RevisedSimplex" can be used for arbitrary- and infinite-precision inputs:

Table[LinearOptimization[obj,N[constr,20],{x,y},Method->met],{met,{"Simplex","RevisedSimplex"}}]

Table[LinearOptimization[obj, constr,{x,y},Method->met],{met,{"Simplex","RevisedSimplex"}}]

Different methods have different strengths, which is typically problem and implementation dependent:

a=SparseArray[{Band[{1,1}]->1.,Band[{1,2}]-> 2.},{500,300}];

"Simplex" and "RevisedSimplex" are good for small dense problems:

Timing[LinearOptimization[Total[x],{a.x\[VectorGreaterEqual]-1,x\[VectorGreaterEqual]0},x,Method->"Simplex"];]

"InteriorPoint" and "CLP" are good for large sparse problems:

Timing[LinearOptimization[Total[x],{a.x\[VectorGreaterEqual]-1,x\[VectorGreaterEqual]0},x,Method->"InteriorPoint"];]

Timing[LinearOptimization[Total[x],{a.x\[VectorGreaterEqual]-1,x\[VectorGreaterEqual]0},x,Method->"CLP"];]

"InteriorPoint" or "CLP" methods may not be able to tell if a problem is infeasible or unbounded:

LinearOptimization[-(x+y),x+y≥1,{x,y}, Method->"InteriorPoint"]

LinearOptimization[-(x+y),x+y≥1,{x,y}, Method->"CLP"]

"Simplex" and "RevisedSimplex" can detect unbounded or infeasible problems:

LinearOptimization[-(x+y),x+y≥1,{x,y}, Method->"Simplex"]

Different methods may give different results for problems with more than one optimal solution:

ContourPlot[-(x+y),{x,y}∈ImplicitRegion[x+y≤1&&x≥0&&y≥0,{x,y}],ContourLabels->True,ImageSize->Small]

"InteriorPoint" may return a solution in the middle of the optimal solution set:

LinearOptimization[-(x+y),{x+y≤1,x≥0,y≥0},{x,y},Method->"InteriorPoint"]

"Simplex" always returns a solution at a corner of the optimal solution set:

LinearOptimization[-(x+y),{x+y≤1,x≥0,y≥0},{x,y},Method->"Simplex"]

OpenersExampleSubsection5(Spacer1)Tolerance  (2)

A smaller (RefLink[Tolerance, paclet:ref/Tolerance]) setting gives a more precise result:

n=2000; m=1000;
(Subscript[x, exact])=ConstantArray[1,n];

Specify the constraint matrices:

a=SparseArray[{Band[{1,1}]->1,Band[{1,2}]-> 2},{m,n}];
b=a.(Subscript[x, exact]);

Find the error between computed and exact minimum value n using different (RefLink[Tolerance, paclet:ref/Tolerance]) settings:

err=Table[{tol,Abs[n-LinearOptimization[Total[x],
{a.x\[VectorGreaterEqual]b,x\[VectorGreaterEqual]1},x,"PrimalMinimumValue",Tolerance->tol,Method->"InteriorPoint"]]},{tol,10.^-Range[7]}]

Visualize the change in minimum value error with respect to tolerance:

ListLogLogPlot[err, PlotRange->All,Joined->True,PlotMarkers->Automatic]

A smaller (RefLink[Tolerance, paclet:ref/Tolerance]) setting gives a more precise answer, but typically takes longer to compute:

n=20000; m=10000;
a=SparseArray[{Band[{1,1}]->1.,Band[{1,2}]-> 2.},{m,n}];
b=Range[m];

A smaller tolerance takes longer:

Timing[ptight=LinearOptimization[Total[x]/n,{a.x\[VectorGreaterEqual]b,x\[VectorGreaterEqual]0},x,"PrimalMinimumValue",Tolerance-> (10^(-12)),Method->"InteriorPoint"];]

Timing[ploose=LinearOptimization[Total[x]/n,{a.x\[VectorGreaterEqual]b,x\[VectorGreaterEqual]0},x,"PrimalMinimumValue",Tolerance-> (10^(-1)),Method->"InteriorPoint"];;]

The tighter tolerance gives a more precise answer:

{ptight,ploose}

OpenersExampleSubsection6(Spacer1)WorkingPrecision  (2)

(RefLink[LinearOptimization, paclet:ref/LinearOptimization]) infers the (RefLink[WorkingPrecision, paclet:ref/WorkingPrecision]) to use from the input:

LinearOptimization[x+y,{x+2y≥3,x≤2,y≤1},{x,y}]

(RefLink[MachinePrecision, paclet:ref/MachinePrecision]) input gives (RefLink[MachinePrecision, paclet:ref/MachinePrecision]) output:

res=LinearOptimization[x+y,{x+2y≥3.,x≤2,y≤1},{x,y}]

Precision[res]

(RefLink[LinearOptimization, paclet:ref/LinearOptimization]) can compute results using higher working precision:

LinearOptimization[x+y, {x+2y≥3,-1≤x≤2,-1≤y≤1},{x,y}, WorkingPrecision-> 30]

If the specified precision is less than the precision of the input arguments, a message is issued:

LinearOptimization[x+y, {x+2y≥3.0,-1≤x≤2,-1≤y≤1},{x,y}, WorkingPrecision->30]

OpenersExampleSection3(Spacer1)Applications  (25)

OpenersExampleSubsection7(Spacer1)Basic Modeling Transformations  (8)

Maximize x+y subject to x+y=3,x≥0,y≤2. Solve a maximization problem by negating the objective function:

LinearOptimization[-(x+y),{x+y==3,x≥0,y≤2},{x,y}]

Negate the primal minimum value to get the corresponding maximal value:

-LinearOptimization[-(x+y),{x+y==3,x≥0,y≤2},{x,y},"PrimalMinimumValue"]

Convert ||a.x-b(Subscript[||, 1]) into a linear function using ||a.x-b(Subscript[||, 1])=(∑_i)(Subscript[s, i]),(Subscript[s, i])=|(Subscript[a, i]).x-(Subscript[b, i])| and -(Subscript[s, i])≤ (Subscript[a, i]).x-(Subscript[b, i])≤ (Subscript[s, i]), i=1,2,..,n:

LinearOptimization[Total[s],{a==, b==Range[10],-s\[VectorLessEqual]a.x-b\[VectorLessEqual]s},{x,s}]

Convert ||a.x-b(Subscript[||, ∞]) into a linear function using ||a.x-b(Subscript[||, ∞])=Max(Subscript[., i])|(Subscript[a, i]).x-(Subscript[b, i])|=s and -s≤ (Subscript[a, i]).x-(Subscript[b, i])≤ s, i=1,2,..,n:

{a,b}=;

LinearOptimization[s,{b-s\[VectorLessEqual]a.x\[VectorLessEqual]b+s},{x,s}]

Minimize |x| subject to the constraints x+y≤3,x≥1,y≥1. Convert |x| into a linear function using |x|=(SuperPlus[x])+(SuperMinus[x]), x=(SuperPlus[x])-(SuperMinus[x]), (SuperPlus[x])≥0,(SuperMinus[x])≥0:

LinearOptimization[xp+xm,{(xp-xm)+y≤3,(xp-xm)≥1,y≥1,xp≥0,xm≥0},{xp,xm,y}]

The problem can also be converted using t=|x|, -t≤x≤t:

LinearOptimization[t,{x+y≤3,x≥1,y≥1,-t≤x≤t},{x,y,t}]

Minimize x+y subject to the constraints |2x+y|≤1,x-y≤3. The constraint can be transformed to -1≤2x+y≤1,x-y≤3:

LinearOptimization[x+y,{-1≤2x+y≤1,x-y≤3},{x,y}]

Minimize x+y subject to the constraints |2x+y|≥1,x-y≤3, x≥0. The constraint |2x+y|≥1 can be interpreted as 2x+y≥1 || 2x+y≤-1. Solve the problem with each constraint:

res1=LinearOptimization[x+y,{2x+y≥1,x-y≤3,x≥0},{x,y},{"PrimalMinimumValue","PrimalMinimizerRules"}]

res2=LinearOptimization[x+y,{2x+y≤-1,x-y≤3,x≥0},{x,y},{"PrimalMinimumValue","PrimalMinimizerRules"}]

The optimal solution is the minimum of the two solutions:

First@MinimalBy[{res1,res2},First]

Show[Plot3D[x+y,],Graphics3D[]]

Minimize f(c.x) subject to a.x⪯b, where f is a nondecreasing function, by instead minimizing c.x. The primal minimizer (SuperStar[x]) will remain the same for both problems. Consider minimizing f(x,y)=(((x+y))^2) subject to 2x-y≤1,x≥ 0,y≥0:

f=Function[{s},s^2];

{minval,(Subscript[X, min])}=LinearOptimization[x+y,{2x-y≥1,x≥0,y≥0},{x,y},{"PrimalMinimumValue","PrimalMinimizerRules"}]

The true minimum value can be obtained by applying the function f to the primal minimum value:

f[minval]

Minimize f(z)=(e^z), z=(1+x+3y), z≥0, subject to 2x+y≥1,-2≤y≤2:

f=Function[{s},Exp[s]];

Since z≥0, solve the problem with 1+x+3y≥0:

{minval,(Subscript[X, min])}=LinearOptimization[(1+x+3y),{2x+y≥1,-2≤y≤2,1+x+3y≥0},{x,y},{"PrimalMinimumValue","PrimalMinimizerRules"}]

The true minimum value can be obtained by applying the function f to the primal minimum value:

f[minval]

OpenersExampleSubsection8(Spacer1)Data-Fitting Problems  (4)

Find a robust linear fit to discrete data by minimizing ||a.x-b(Subscript[||, 1])=(∑_i)|(Subscript[a, i]).x-(Subscript[b, i])|. Using auxiliary variables (Subscript[s, i])=|(Subscript[a, i]).x-(Subscript[b, i])|, the problem is transformed to minimize (∑_i)(Subscript[s, i]) subject to the constraints -(Subscript[s, i])≤ (Subscript[a, i]).x-(Subscript[b, i])≤ (Subscript[s, i]), which is a linear optimization problem:

{input,output}=;
dataplot=ListPlot[output,PlotRange->All]

Fit the data:

res=LinearOptimization[Total[s],{a==input,b==output,-s\[VectorLessEqual]a.x-b\[VectorLessEqual]s},{x,s}]

A plot of the line shows that it is robust to the outliers:

Show[ dataplot,Plot[Evaluate[x.{1,t}/.res],{t,1,30}]]

Compare with the fit obtained using (RefLink[LeastSquares, paclet:ref/LeastSquares]):

line = Fit[output,{1,t},t]

Show[ dataplot,Plot[Evaluate[{x.{1,t}/.res,line}],{t,1,30}]]

Find a linear fit by minimizing ||a.x-b(Subscript[||, 1]) such that the line should pass through the first data point:

{input,output}=;

By adding the constraint a[[1]].x==b[[1]], the line is forced to pass through the first point:

res=LinearOptimization[Total[s],{a==input,b==output,-s\[VectorLessEqual]a.x-b\[VectorLessEqual]s,input[[1]].x==output[[1]]},{x,s}]

The plot shows that the line passes through the first data point while still maintaining a robust fit:

Show[Plot[Evaluate[x.{1,t}/.res],{t,1,30}, Epilog->{Red,PointSize[0.02],Point[{1,output[[1]]}]}], ListPlot[output], PlotRange->All,AxesOrigin-> {0,0}]

Find a robust fit to nonlinear discrete data by minimizing ||a.x-b(Subscript[||, 1]). Generate some noisy data:

data=;
dataPlot = ListPlot[data, PlotRange->All]

Fit the data using the bases (Subscript[ϕ, i])(x)=(((1+(((x-(Subscript[c, i])))^2)))^(1/2)). The approximating function will be s(x)=(∑_(i=1))(Subscript[λ, i])(Subscript[ϕ, i])(x):

bases = Sqrt[1+(x-#)^2]&/@Subdivide[-3,3,10];
input=DesignMatrix[data,bases,x,IncludeConstantBasis->False];
output=data[[All,2]];

Find the coefficients of the approximating function:

res=LinearOptimization[Total[s],{a==input,b==output,-s\[VectorLessEqual]a.λ-b\[VectorLessEqual]s}, {λ,s}];

Visualize the fit:

Show[dataPlot,Plot[Evaluate[bases.λ/.res],{x,-3,3}]]

Compare the approximating model with the reference and the least squares fit:

lsqfit=Fit[data,bases,x];

Plot[{Sin[2x]+Cos[x^2],Evaluate[bases.λ/.res],lsqfit},{x,-3,3}, PlotLegends->{"Reference","(Subscript[L, 1])fit","Least squares fit"}]

Find a linear fit to discrete data by minimizing ||a.x-b(Subscript[||, ∞]). Since ||a.x-b(Subscript[||, ∞])=Max(Subscript[., i])|(Subscript[a, i]).x-(Subscript[b, i])|, by introducing the scalar auxiliary variable s=Max.|(Subscript[a, i]).x-(Subscript[b, i])|, the problem can be transformed to minimize s, subject to the constraints -s≤ (Subscript[a, i]).x-(Subscript[b, i])≤ s, i=1,2,..,n:

{a,b}=;
dataplot=ListPlot[b,PlotRange->All]

Fit the data:

res = LinearOptimization[s,{-s+b\[VectorLessEqual]a.x\[VectorLessEqual]s+b},{x,s}]

The line represents the "best-worst case" fit:

Show[dataplot,Plot[x.{1,t}/.res,{t,1,30}]]

The function (RefLink[Fit, paclet:ref/Fit]) can be used directly to compute the coefficients:

Fit[b,{1,t},t,NormFunction->Infinity]

OpenersExampleSubsection9(Spacer1)Classification Problems  (3)

Find a line a.x+b that separates two groups of points P={(Subscript[p, 1]),(Subscript[p, 2]),..(Subscript[p, n])} and Q={(Subscript[q, 1]),(Subscript[q, 2]),..(Subscript[q, m])}:

{set1,set2}=;
dataPlot = ListPlot[{set1,set2},PlotMarkers->{Automatic,PointSize[0.025]}]

For separation, set 1 must satisfy a.(Subscript[p, i])+b \[VectorGreaterEqual]1 and set 2 must satisfy a.(Subscript[q, j])+b\[VectorLessEqual]-1:

res=LinearOptimization[0,{set1.a+b\[VectorGreaterEqual]1,set2.a+b\[VectorLessEqual]-1},{a,b}]

The separating line is:

line=a.{x,y}+b /. res

Show[RegionPlot[line≤0,{x,0,1},{y,0,1}],dataPlot]

Find a quadratic polynomial (Subscript[c, 0])+(Subscript[c, 1])x+(Subscript[c, 2])y+(Subscript[c, 3])(x^2)+(Subscript[c, 4])(y^2)+(Subscript[c, 5])x y that separates two groups of points P={(Subscript[p, 1]),(Subscript[p, 2]),..(Subscript[p, n])} and Q={(Subscript[q, 1]),(Subscript[q, 2]),..(Subscript[q, m])}:

{set1,set2}=;
dataPlot = ListPlot[{set1,set2},PlotMarkers->{Automatic,PointSize[0.025]}]

For separation, set 1 must satisfy (Subscript[c, 0])+(Subscript[c, 1])(Subscript[p, x,i])+(Subscript[c, 2])(Subscript[p, y,i])+(Subscript[c, 3])((Subscript[p, x,i])^2)+(Subscript[c, 4])((Subscript[p, y,i])^2)+(Subscript[c, 5])(Subscript[p, x,i]) (Subscript[p, y,i])≥1 and set 2 must satisfy (Subscript[c, 0])+(Subscript[c, 1])(Subscript[q, x,j])+(Subscript[c, 2])(Subscript[q, y,j])+(Subscript[c, 3])((Subscript[q, x,j])^2)+(Subscript[c, 4])((Subscript[q, y,j])^2)+(Subscript[c, 5])(Subscript[q, x,j]) (Subscript[q, y,j])≤1:

poly[{x_,y_}]:={1,x,y,(x^2),(y^2),x y};
set1Quad=poly/@set1;
set2Quad = poly/@set2;

The coefficients of the quadratic polynomial are:

res=LinearOptimization[0,{set1Quad.c\[VectorGreaterEqual]1,set2Quad.c\[VectorLessEqual]-1},c]

The quadratic polynomial is:

quad=poly[{x,y}].(c/. res)

Show[RegionPlot[quad≤0,{x,-1/2,3/2},{y,-1/2,3/2}],dataPlot]

Find the minimum-degree polynomial that can separate two sets of points in the plane:

set1=N[{{0,0},{-2,-2},{2,-2},{2,2},{-2,2}}];
set2=N[{{0,-1},{1,0},{0,1},{-1,0}}];

Define a polynomial power function that avoids issues when x or y is 0:

pPower[{x_,y_},{i_,j_}]:=If[i==0,1,(x^i)]If[j==0,1,(y^j)];

Define a function that gives a polynomial of degree d with coefficient a as a function of {x,y}:

poly[d_][{x_,y_}]:=Sum[(Subscript[a, i,j])pPower[{x,y},{i,j}],{j,0,d},{i,0,d-j}];

The variables for a degree d are the coefficients:

variables[d_]:=Flatten[Table[(Subscript[a, i,j]),{j,0,d},{i,0,d-j}]];

The constraints enforce separation between set 1 and set 2:

constraints[d_]:=Join[Thread[Map[poly[d],set1]≤-1], Thread[Map[poly[d],set2]≥1]];

For example, the constraints for quadratics are:

TableForm[constraints[2]]

For separation, the only condition is that all the constraints must be satisfied. The objective vector is set to 0:

d=0;
maxDegree=10;
minimizer=Indeterminate;
While[Not[VectorQ[minimizer,NumberQ]],
If[++d>maxDegree,d=$Failed; Break[]];
minimizer= Quiet[LinearOptimization[0,constraints[d],
variables[d], "PrimalMinimizer"]]
];
d

Find a separating polynomial and show it with the points:

rules = LinearOptimization[0,constraints[d],variables[d], "PrimalMinimizerRules"];
separator = poly[d][{x,y}] /. rules

RegionPlot[ImplicitRegion[Evaluate[separator≤ 0],{x,y}],Epilog->{PointSize[0.025],Point[set1],Red,Point[set2]},PlotRange->{{-3,3},{-3,3}}]

OpenersExampleSubsection10(Spacer1)Manufacturing Problems  (2)

Find the mix of products to manufacture that maximizes the profits of a company. The company makes three products. The revenue per unit, cost per unit and maximum capacity are given by:

revenuePerUnit={90.,100.,70.};
costPerUnit={45.,40.,20.};
capacity= {100,40,60};

Each product is made using a combination of four machines. The time a machine spends on each product is:

machineTime ={{20,10,10},{12,28,16},{15,6,5},{10,15,0}};

Row[BarChart[machineTime[[#]],]&/@Range[4]]

The profit equals revenue minus cost times the number of units (Subscript[x, i]) for product i:

profit = (revenuePerUnit - costPerUnit ).x

Each machine can run at most 2400 minutes per week:

timeConstraint=machineTime.x\[VectorLessEqual]2400;

The company makes between zero and the capacity for each product:

capacityConstraint=0\[VectorLessEqual]x\[VectorLessEqual]capacity;

The optimal product mix is given by:

res=LinearOptimization[-profit,{timeConstraint,capacityConstraint},x]

The resulting profit is given by:

profit/.res

The breakdown of machining time for each product is:

GraphicsRow[, ImageSize->450]

The same company wants find out what capacity to change that will have the largest impact on profit. The optimal product mix and profit as formulated previously are given by:

timeConstraint={{20,10,10},{12,28,16},{15,6,5},{10,15,0}}.x\[VectorLessEqual]2400;
capacity=x\[VectorLessEqual]{100,40,60};
profit = {45.,60.,50.}.x;

The profit based on the current capacity is:

max0 = -LinearOptimization[-profit,{timeConstraint,capacity},x,"PrimalMinimumValue"]

Extract all the data needed to analyze the impact of the capacity change on profit:

{sens,{a,b}}=LinearOptimization[-profit,{timeConstraint,capacity} ,x,
{"ConstraintSensitivity","LinearInequalityConstraints"}]

The sensitivity gives the derivative of the profit wrt each of the constraints. Extract the constraints that have any sensitivity at all:

sensPos=Flatten@Position[Chop[sens[[1]]],_?(Abs[#]>0&)]

Extract the sensitive constraints from a.x+b⪰0:

Table[a[[i]].{(Subscript[x, 1]),(Subscript[x, 2]),(Subscript[x, 3])}+b[[i]]≥0,{i,sensPos}]//TableForm

If the production capacity for product 3 is increased, there should be an increase in profit:

capacity1=x\[VectorLessEqual]{100,400,60+10};

max1=-LinearOptimization[-profit,{timeConstraint,capacity1} ,x,"PrimalMinimumValue"]

The resulting increase in profit:

max1-max0

This corresponds exactly to the predicted increase from the sensitivity:

sens[[1,7]]×10

If the production capacity of product 1, product 2 or both is changed, the profit does not increase:

capacity2=x\[VectorLessEqual]{100+10,40+10,60};

max2=-LinearOptimization[-profit,{timeConstraint,capacity2},x,"PrimalMinimumValue"];
max2-max0

OpenersExampleSubsection11(Spacer1)Transportation Problems  (1)

John Conner has learned that Skynet is planning on attacking five settlements. He needs to deploy soldiers from three base camps to the settlements to meet the threat. The fuel cost to transport one soldier from base camp i to settlement (Subscript[s, j]) is:

fuelCost={{10,15,15,10,20},{15,10,10,15,15},{10,20,5,10,10}};

GraphicsRow[BarChart[fuelCost[[#]],]&/@Range[3], ]

Let (Subscript[r, i,j]) represent the number of soldiers to be shipped from base camp i to settlement j. The cost function to minimize is:

totalCost=Total[fuelCost*r,2];

Each base camp needs to distribute a maximum of 100 soldiers between the settlements:

maxSoldiersDeployed=Total[r,{2}]\[VectorLessEqual]100;

Each settlement needs at least 50 soldiers from a combination of the base camps:

minSoldiersNeeded={Total[r]\[VectorGreaterEqual]50,r\[VectorGreaterEqual]0};

The number of soldiers to be deployed from each base camp to each settlement is:

res=LinearOptimization[totalCost,{maxSoldiersDeployed,minSoldiersNeeded},r∈Matrices[Dimensions[fuelCost]]]

The breakdown of soldier deployment from each base camp is:

TableForm[r/.res, TableHeadings->{{"Base Camp 1","Base camp 2","Base camp 3"},Array[(Subscript[s, #])&,5]}]

OpenersExampleSubsection12(Spacer1)Optimal Assignment Problems  (2)

The Jedi Order is fighting the separatist armies led by Dooku. Dooku and two of his generals, Asajj Ventress and General Grievous, have attacked three outer colonies. The Jedi Order sends three Jedi to assist the outer colonies. Let (Subscript[p, ij]) be the odds of Jedi i beating Sith j:

p={{0.4,0.5,0.5},{0.6,0.8,0.6},{0.7,0.9,0.6}};

GraphicsRow[BarChart[p[[#]],]&/@Range[3], ]

Let x be a matrix variable where (Subscript[x, ij]) is 1 if Jedi i battles Sith j, and 0 otherwise. The cost function is the expected number of defeats given by (∑_(i=1))(∑_(j=1))(1-(Subscript[p, ij]))(Subscript[x, ij]) and constructed using (RefLink[Inactive, paclet:ref/Inactive])[(RefLink[Times, paclet:ref/Times])]:

defeats=Total[(1-p)*x,2];

Since the outer colonies are far apart, each Jedi may only be assigned to one Sith, and vice versa. Also, since (Subscript[x, ij]) is 0 or 1, then (Subscript[x, ij])≥ 0:

constraints={Total[x] - 1== 0, Total[x,{2}] -1==0, x\[VectorGreaterEqual]0};

Find the appropriate Jedi to face the Sith to minimize the chances of defeat:

rule=LinearOptimization[defeats, constraints,x∈Matrices[Dimensions[p]]]

This means that Obi-Wan Kenobi is sent to General Grievous, Mace Windu is sent to Asajj Ventress and Ki-Adi-Mundi is sent to Dooku. The expected number of defeats is:

LinearOptimization[defeats, constraints,x∈Matrices[Dimensions[p]],"PrimalMinimumValue"]

Solve a bipartite matching problem between two sets a and b that minimizes the total weight:

weights=;

Let x be a matrix variable, where (Subscript[x, ij]) is 1 when (Subscript[a, i]) is matched to (Subscript[b, j]) and 0 otherwise. The total weight of the matching is given by (∑_(i=1))(∑_(j=1))(Subscript[w, ij])(Subscript[x, ij]) and constructed using (RefLink[Inactive, paclet:ref/Inactive])[(RefLink[Times, paclet:ref/Times])]:

cost=Total[weights*x,2];

The row and column sums of x need to be 1 and also x≥0:

constraints={Total[x]==1,Total[x,{2}]==1,x\[VectorGreaterEqual]0};

Use (RefLink[LinearOptimization, paclet:ref/LinearOptimization]) to solve the problem:

{minCost, rule} = LinearOptimization[cost, constraints, x∈Matrices[Dimensions[weights]],{"PrimalMinimumValue", "PrimalMinimizerRules"}]

The matching permutation may be extracted by converting the matrix x into a graph. (RefLink[Round, paclet:ref/Round]) is used to truncate small values that may occur due to numerical error in some methods:

graph=AdjacencyGraph[Round[x/.rule]];

Use the adjacency graph to build a bipartite graph:

edges=EdgeList[graph] /. DirectedEdge[i_,j_]:> DirectedEdge[(Subscript[a, i]),(Subscript[b, j])];
Graph[edges,]

OpenersExampleSubsection13(Spacer1)Inventory Control Problems  (2)

Find the inventory that a retail store must order for a product P per week to minimize cost. Let (Subscript[x, i]) be the stock available at the beginning of the week i. Let (Subscript[u, i]) be the stock order received from the wholesaler at the beginning of the week i. The cost to purchase one unit from the wholesaler is $3:

purchaseCost=3u;
purchaseConstraint = u\[VectorGreaterEqual]0;

The cost of holding excess or negative inventory is $1 per unit. The inventory cost is |x|, which is not linear but can be represented as |(Subscript[x, i])|=(Subscript[((SuperPlus[x])), i])+(Subscript[((SuperMinus[x])), i]), (Subscript[x, i])=(Subscript[((SuperPlus[x])), i])-(Subscript[((SuperMinus[x])), i]), (Subscript[((SuperPlus[x])), i])≥0 ,(Subscript[((SuperMinus[x])), i])≥0:

inventoryCost=xp+xm;
inventoryConstraints={xp\[VectorGreaterEqual]0,xm\[VectorGreaterEqual]0};

The total cost to minimize is (∑_(i=1))|(Subscript[x, i])| +3(Subscript[u, i]):

totalCost=Total[inventoryCost+purchaseCost];

Let (Subscript[w, i]) be the demand at the beginning of the week i. The demand over a period of 52 weeks is given by 5 Sin(2 π t/52)+6 for 0≤t≤52. The store inventory evolves as (Subscript[x, i+1])=(Subscript[x, i])+(Subscript[u, i])-(Subscript[w, i]), i=0,1,2..,52:

n=52;
w=6.+5Sin[Subdivide[0,2Pi,n]];
a=SparseArray[{Band[{1,1}]->-1., Band[{1,2}]->1.},{n+1,n+1}];
demandConstraint = a.(xp-xm) == Inactive[Plus][u,-w];

At week 0, there is an existing inventory (Subscript[x, 0])=10 units and stock order (Subscript[u, 0])=5 units. The store cannot purchase more than 10 units per week from the wholesaler:

initialInventory={Indexed[xp,1]-Indexed[xm,1]==10, Indexed[u,1]==5};
maxOrder=u\[VectorLessEqual]10;

Find the optimal inventory to keep in store and order from the wholesaler by minimizing cost:

vars= Element[#,Vectors[n+1,Reals]]&/@{xp,xm,u};
res=LinearOptimization[totalCost,{purchaseConstraint,inventoryConstraints,demandConstraint,initialInventory,maxOrder},vars];

Here is a plot of inventory, stock order and demand that shows that after week 26, the store does not need to carry any inventory, hence minimizing inventory cost:

ListLinePlot[{xp-xm,u, w} /. res, PlotRange->All, PlotStyle->{Red,Blue,Green}, PlotLegends->{"Inventory","Stock order","Demand"}]

Find the inventory that a retail store must order for a product P per week to minimize cost without back orders. Let (Subscript[u, i]) be the stock order received from the wholesaler at the beginning of the week i. The cost of purchasing one unit from the wholesaler is $3:

purchaseCost=3u;
purchaseConstraint = u\[VectorGreaterEqual]0;

The cost of holding excess or negative inventory is $1 per unit:

inventoryCost=xp+xm;
inventoryConstraints={xp\[VectorGreaterEqual]0,xm\[VectorGreaterEqual]0};

Let (Subscript[w, i]) be the demand at the beginning of the week i. The demand over a period of 52 weeks is given by 5 Sin(2 π t/52)+6 for 0≤t≤52. The store inventory evolves as (Subscript[x, i+1])=(Subscript[x, i])+(Subscript[u, i])-(Subscript[w, i]), i=0,1,2..,52:

n=52;
w=6.+5Sin[Subdivide[0,2Pi,n]];
a=SparseArray[{Band[{1,1}]->-1., Band[{1,2}]->1.},{n+1,n+1}];
demandConstraint = a.(xp-xm) == Inactive[Plus][u,-w];

The initial inventory is (Subscript[x, 0])=10  units and (Subscript[u, 0])=5 units, and the maximum purchase from the wholesaler is 10 units:

initialInventory={Indexed[xp,1]-Indexed[xm,1]==10, Indexed[u,1]==5};
maxOrder=u\[VectorGreaterEqual]0;

The company wants to make sure that there are no back orders, i.e. (Subscript[x, i])≥0:

noBackOrders=xp-xm\[VectorGreaterEqual]0;

Find the optimal inventory to keep in store and order from the wholesaler by minimizing cost:

vars= Element[#,Vectors[n+1,Reals]]&/@{xp,xm,u};
res = LinearOptimization[Total[purchaseCost+inventoryCost], {purchaseConstraint, inventoryConstraints,demandConstraint,initialInventory,maxOrder,noBackOrders}, vars];

This plot shows that in order to avoid back orders, excess inventory must be purchased earlier in the weeks to meet the demand. The store has no inventory costs by weeks 20–21:

ListLinePlot[{xp-xm,u, w} /. res, PlotRange->All, PlotStyle->{Red,Blue,Green}, PlotLegends->{"Inventory","Stock order","Demand"}]

OpenersExampleSubsection14(Spacer1)Structural Optimization Problems  (2)

Minimize the weight of a two-beam truss subject to stress and displacement constraints:

The weight is W=ρ ((Subscript[A, 1])(Subscript[l, 1])+(Subscript[A, 2])(Subscript[l, 2])), where ρ is the density, and (Subscript[A, 1]) and (Subscript[A, 2]) are the cross-sectional areas of bar 1 and bar 2. For this problem, assume (Subscript[A, 1])=(Subscript[A, 2])=A:

weight=ρ A l(1/Cos[α]+1);

The forces acting along the axis of bar 1 and bar 2 are:

{(Subscript[F, 1]),(Subscript[F, 2])}={F/Sin[α],F/Tan[α]};

The stress σ=(Subscript[F, i])/(Subscript[A, i]), where (Subscript[F, i]) is the force acting along the bar, must be less than a specified stress (Subscript[σ, max]):

forceConstraints={(Subscript[F, 1])≤A (Subscript[σ, max]),(Subscript[F, 2])≤A (Subscript[σ, max])};

The displacement at the free node is ((Subscript[δ, x]),(Subscript[δ, y])). The component of this displacement along the axis of the bar must equal the displacement (Subscript[δ, i])=((Subscript[F, i])(Subscript[l, i]))/(Y (Subscript[A, i])), where Y is the Young's modulus:

deqns={((Subscript[F, 1])l/Cos[α])/(Y A)==(Subscript[δ, x])/Cos[α]+(Subscript[δ, y])/Sin[α],((Subscript[F, 2]) l)/(Y A)==(Subscript[δ, x])};

The vertical displacement of the free node must be less than a given δ:

displacementConstraint=A((Subscript[δ, y])/.First@Solve[deqns,{(Subscript[δ, x]),(Subscript[δ, y])}])≤δ A;

The material being used is aluminium. The design parameters are:

params={Y->69 (10^6),ρ->2700.,(Subscript[σ, max])->69000,α->Pi/6,l->1,F->(10^3),δ->(10^(-2))};

Compute the minimum weight and the optimal cross-sectional areas of the beam:

LinearOptimization[weight/.params,{forceConstraints, displacementConstraint}/.params,A,{"PrimalMinimumValue","PrimalMinimizerRules"}]

Find the range of angle and length for which the truss weight does not exceed 180 kg:

minTrussWeight[α_?NumericQ,l_?NumericQ]:=Module[{A},LinearOptimization[,A,"PrimalMinimumValue"]];

ContourPlot[minTrussWeight[x,y],]

Minimize the weight of a two-beam truss subject to stress and displacement constraints, and find the constraints that have a major effect on the design of the truss:

The weight is W=ρ ((Subscript[A, 1])(Subscript[l, 1])+(Subscript[A, 2])(Subscript[l, 2])), where ρ is the density, and (Subscript[A, 1]) and (Subscript[A, 2]) are the cross-sectional areas of bar 1 and bar 2. For this problem, assume (Subscript[A, 2])=β (Subscript[A, 1]):

weight=ρ A l(1/Cos[α]+β);

The constraint due to the forces acting on the bars is:

{(Subscript[F, 1]),(Subscript[F, 2])}={F/Sin[α],F/Tan[α]};
forceConstraints={A (Subscript[σ, max])-(Subscript[F, 1])≥ 0,A β (Subscript[σ, max])-(Subscript[F, 2])≥ 0};
displacementConstraint=δ A β Y-F l(β Sec[α]-1)≥0;

The material being used is aluminium. The design parameters are:

params={Y->69 (10^6),ρ->2700.,(Subscript[σ, max])->69000,α->Pi/6,l->1,F->(10^3),δ->(10^(-2)),β->1.1};

Minimize the weight of the truss by finding the optimal cross-sectional area A of the bars:

LinearOptimization[weight/.params,{forceConstraints, displacementConstraint}/.params,A,{"PrimalMinimumValue","PrimalMinimizerRules"}]

Large constraint sensitivity values indicate a greater influence on the design:

{sens,{a,b}}=LinearOptimization[weight/.params,{forceConstraints, displacementConstraint}/.params,A,{"ConstraintSensitivity", "LinearInequalityConstraints"}];
sens

si=First@Flatten@Position[Chop[sens[[1]]],x_/;Abs[x]>0]

The critical constraint influencing the design is the force constraint (Subscript[F, 1])/A≤σ acting on bar 1:

{a[[si]].{A}+b[[si]] ≥ 0, forceConstraints[[1]] /.params}

OpenersExampleSubsection15(Spacer1)Sequential Linear Optimization  (1)

Minimize a nonlinear function f(x) subject to nonlinear constraints g(x)≥0. The minimization can be done by approximating f (x) as f(x+d) ≈f(x)+ (∇f).d and the constraints as g(x) as g(x+d)=g(x)+(∇g).d≥0. This leads to a linear minimization problem that can be solved iteratively. Consider a case where f=-(x^2)Exp(-(y^2)) and g={(x^2)+(y^2)≤25, x+2y-12≤0,x≥0,y≥0}:

f[{x_,y_}]:=-(x^2) Exp[-(y^2)];
g[{x_,y_}]:={25-(x^2)-(y^2),12-x-2y,x,y};

The gradients of the minimizing function and constraints are:

df[{x_,y_}]=D[f[{x,y}],{{x,y}}];
dg[{x_,y_}]=D[g[{x,y}],{{x,y}}];

The subproblem is to find d by minimizing (∇f((SuperStar[x]))).d subject to (Subscript[g, i])((SuperStar[x]))+(∇(Subscript[g, i])((SuperStar[x]))). d≥0:

LinearSubProblem[xiter_]:= Module[{c,a,b},
{c,a,b}=#[xiter]&/@{df,dg,g};
LinearOptimization[c,{a,b}]
];

Iterate, starting with an initial guess of (SuperStar[x])=(2, 1). The next iterate is (Subscript[x, new])=(SuperStar[x])+s.d, where s∈[0,1] is the step length to ensure that the constraints are always satisfied:

(SuperStar[x])={2.,1.};
xp = Last@Reap[Do[
Sow[(SuperStar[x]),"IterationPoints"];
d = LinearSubProblem[(SuperStar[x])];
If[Norm[d]≤10^-5, Break[]];
s=1;
While[Min[g[(SuperStar[x])+s*d]]≤0,s/=2];
(SuperStar[x])+=s*d;
,{20}],"IterationPoints"]

Visualize the convergence of the result. The final result is the green point:

ContourPlot[f[{x,y}],{x,0,6},{y,0,2},]

OpenersExampleSection4(Spacer1)Properties & Relations  (9)

(RefLink[LinearOptimization, paclet:ref/LinearOptimization]) gives the global minimum of the objective function:

res=LinearOptimization[x+2y,{3x+4y≥5,x≥0,y≥0},{x,y}]

Show[Plot3D[x+2y,{x,y}∈ImplicitRegion[3x+4y≥5&&x≥0&&y≥0,{x,y}]],Graphics3D[{PointSize[0.05],Point[{x,y,x+2y}/.res]}]]

(RefLink[Minimize, paclet:ref/Minimize]) also gives global exact results for linear optimization problems:

Minimize[{x+2y,{3x+4y≥5&&x≥0&&y≥0}},{x,y}]

(RefLink[NMinimize, paclet:ref/NMinimize]) can be used to obtain approximate results using global methods:

NMinimize[{x+2y,{3x+4y≥5&&x≥0&&y≥0}},{x,y}]

(RefLink[FindMinimum, paclet:ref/FindMinimum]) can be used to obtain approximate results using local methods:

FindMinimum[{x+2y,{3x+4y≥5&&x≥0&&y≥0}},{x,y}]

(RefLink[LinearFractionalOptimization, paclet:ref/LinearFractionalOptimization]) uses (RefLink[LinearOptimization, paclet:ref/LinearOptimization]) for solving problems:

LinearFractionalOptimization[((Subscript[x, 1])+(Subscript[x, 2]))/(1+3(Subscript[x, 1])+2(Subscript[x, 2])),{(Subscript[x, 1])+(Subscript[x, 2])≥2,(Subscript[x, 1])≥0,(Subscript[x, 2])≥0},{(Subscript[x, 1]),(Subscript[x, 2])}]

Transform the objective and constraints using t=1/(1+3(Subscript[x, 1])+2(Subscript[x, 2])),(Subscript[y, 1])=(Subscript[x, 1])t,(Subscript[y, 2])=(Subscript[x, 2])t:

obj=((Subscript[y, 1])+(Subscript[y, 2]));
constraints={(Subscript[y, 1])+(Subscript[y, 2])-2t≥0,(Subscript[y, 1])≥0,(Subscript[y, 2])≥0,t≥0,3(Subscript[y, 1])+2(Subscript[y, 2])+t==1};

res=LinearOptimization[obj,constraints,{(Subscript[y, 1]),(Subscript[y, 2]),t}]

Obtain the solution to the original problem by applying the inverse transformation:

{(Subscript[y, 1])/t,(Subscript[y, 2])/t}/.res

(RefLink[QuadraticOptimization, paclet:ref/QuadraticOptimization]) is a generalization of (RefLink[LinearOptimization, paclet:ref/LinearOptimization]):

LinearOptimization[x+2y,{3x+4y≥5,x≥0,y≥0},{x,y}]

QuadraticOptimization[x+2y,{3x+4y≥5,x≥0,y≥0},{x,y}]

(RefLink[SecondOrderConeOptimization, paclet:ref/SecondOrderConeOptimization]) is a generalization of (RefLink[LinearOptimization, paclet:ref/LinearOptimization]):

LinearOptimization[x+2y,{3x+4y≥5,x≥0,y≥0},{x,y}]

SecondOrderConeOptimization[x+2y,{3x+4y≥5,x≥0,y≥0},{x,y}]

(RefLink[SemidefiniteOptimization, paclet:ref/SemidefiniteOptimization]) is a generalization of (RefLink[LinearOptimization, paclet:ref/LinearOptimization]):

LinearOptimization[x+2y,{3x+4y≥5,x≥0,y≥0},{x,y}]

SemidefiniteOptimization[x+2y,{3x+4y≥5,x≥0,y≥0},{x,y}]

(RefLink[ConicOptimization, paclet:ref/ConicOptimization]) is a generalization of (RefLink[LinearOptimization, paclet:ref/LinearOptimization]):

LinearOptimization[x+2y,{3x+4y≥5,x≥0,y≥0},{x,y}]

ConicOptimization[x+2y,{3x+4y≥5,x≥0,y≥0},{x,y}]

OpenersExampleSection5(Spacer1)Possible Issues  (3)

Constraints specified using strict inequalities may not be satisfied:

res=LinearOptimization[x+y,{x-2y≤3,x>1,y>0},{x,y}]

The reason is that (RefLink[LinearOptimization, paclet:ref/LinearOptimization]) solves inf {x+y|x-2 y≤3,x>1,y>0}:

{x-2y≤3,x>1,y>0}/.res

The minimum value of an empty set or infeasible problem is defined to be ∞:

LinearOptimization[x,{x≤0,x≥1},x,"PrimalMinimumValue"]

The minimizer is (RefLink[Indeterminate, paclet:ref/Indeterminate]):

LinearOptimization[x,{x≤0,x≥1},x]

The minimum value for an unbounded set or unbounded problem is -∞:

LinearOptimization[-(x+y),x+y≥1,{x,y},"PrimalMinimumValue"]

The minimizer is (RefLink[Indeterminate, paclet:ref/Indeterminate]):

LinearOptimization[-(x+y),x+y≥1,{x,y}]

({{DocumentationFooterIconsRelatedFunction.png, ({{See Also}, {(RefLink[LinearFractionalOptimization, paclet:ref/LinearFractionalOptimization])  ▪  (RefLink[QuadraticOptimization, paclet:ref/QuadraticOptimization])  ▪  (RefLink[SecondOrderConeOptimization, paclet:ref/SecondOrderConeOptimization])  ▪  (RefLink[SemidefiniteOptimization, paclet:ref/SemidefiniteOptimization])  ▪  (RefLink[ConicOptimization, paclet:ref/ConicOptimization])}})}})

({{DocumentationFooterIconsRelatedGuide.png, ({{Related Guides}, {▪(RefLinkPlain[Optimization, paclet:guide/Optimization])}, {▪(RefLinkPlain[Convex Optimization, paclet:guide/ConvexOptimization])}})}})