Nonlinear Programming¶
Problems that fit the general linear programming format but contain nonlinear functions are termed nonlinear programming problems.
Solution methods are more complex than linear programming methods.
Determining an optimal solution is often difficult, if not impossible.
Solution techniques generally involve searching a solution surface for high or low points requiring the use of advanced mathematics.
A nonlinear problem containing one or more constraints becomes a constrained optimization model or a nonlinear programming model.
A nonlinear programming model has the same general form as the linear programming model except that the objective function and/or the constraint(s) are nonlinear.
Solution procedures are much more complex and no guaranteed procedure exists for all nonlinear problem models.
Unlike linear programming, the solution is often not on the boundary of the feasible solution space.
Cannot simply look at points on the solution space boundary but must consider other points on the surface of the objective function.
This greatly complicates solution approaches.
Solution techniques can be very complex.
Nonlinear problem example:¶
Solve the following nonlinear programming model:
\[\begin{split}
\begin{align}
&\text{max}\\
&\qquad z=(4-0.1x_1)x_1+(5-0.2x_2)x_2\\
&\text{s.t.}\\
&\qquad x_1+2x_2=40\\
&\qquad x_i \ge 0
\end{align}
\end{split}\]
!pip install gurobipy
# Import gurobi library
from gurobipy import * # This command imports the Gurobi functions and classes.
import numpy as np
# create new model
m = Model('Beaver Creek Pottery Company ')
x1 = m.addVar(lb = 0, vtype = GRB.CONTINUOUS, name='x1')
x2 = m.addVar(lb = 0, vtype = GRB.CONTINUOUS, name='x2')
m.setObjective((4-0.1*x1)*x1+(5-0.2*x2)*x2, GRB.MAXIMIZE)
m.addConstr(1*x1+2*x2==40, 'const1')
m.optimize()
#display optimal production plan
for v in m.getVars():
print(v.varName, v.x)
print('optimal total revenue:', m.objVal)
Collecting gurobipy
Downloading gurobipy-9.5.0-cp37-cp37m-manylinux2014_x86_64.whl (11.5 MB)
|████████████████████████████████| 11.5 MB 4.7 MB/s
?25hInstalling collected packages: gurobipy
Successfully installed gurobipy-9.5.0
Restricted license - for non-production use only - expires 2023-10-25
Gurobi Optimizer version 9.5.0 build v9.5.0rc5 (linux64)
Thread count: 1 physical cores, 2 logical processors, using up to 2 threads
Optimize a model with 1 rows, 2 columns and 2 nonzeros
Model fingerprint: 0x34cdafe3
Model has 2 quadratic objective terms
Coefficient statistics:
Matrix range [1e+00, 2e+00]
Objective range [4e+00, 5e+00]
QObjective range [2e-01, 4e-01]
Bounds range [0e+00, 0e+00]
RHS range [4e+01, 4e+01]
Presolve time: 0.02s
Presolved: 1 rows, 2 columns, 2 nonzeros
Presolved model has 2 quadratic objective terms
Ordering time: 0.00s
Barrier statistics:
AA' NZ : 0.000e+00
Factor NZ : 1.000e+00
Factor Ops : 1.000e+00 (less than 1 second per iteration)
Threads : 1
Objective Residual
Iter Primal Dual Primal Dual Compl Time
0 -5.98930000e+05 6.12060000e+05 2.00e+03 0.00e+00 1.01e+06 0s
1 7.00620769e+01 1.96612266e+04 2.00e-03 0.00e+00 9.80e+03 0s
2 7.00738420e+01 9.86563122e+01 9.19e-07 0.00e+00 1.43e+01 0s
3 7.04164770e+01 7.06320152e+01 2.71e-10 0.00e+00 1.08e-01 0s
4 7.04166667e+01 7.04168819e+01 0.00e+00 8.88e-16 1.08e-04 0s
5 7.04166667e+01 7.04166669e+01 0.00e+00 0.00e+00 1.08e-07 0s
6 7.04166667e+01 7.04166667e+01 0.00e+00 8.88e-16 1.08e-10 0s
Barrier solved model in 6 iterations and 0.09 seconds (0.00 work units)
Optimal objective 7.04166667e+01
x1 18.33333333336888
x2 10.833333333315561
optimal total revenue: 70.41666666666666
A Nonlinear Programming Model with Multiple Constraints:¶
Solve the following nonlinear programming model:
\[\begin{split}
\begin{align}
&\text{max}\\
&\qquad z=\left(\frac{1500-x_1}{24.6}-12\right)x_1
+\left(\frac{2700-x_2}{63.8}-9\right)x_2\\
&\text{s.t.}\\
&\qquad 2x_1+2.7x_2\le 6000\\
&\qquad 3.6x_1+2.9x_2\le 8500\\
&\qquad 7.2x_1+8.5x_2\le 15000\\
&\qquad x_i \ge 0
\end{align}
\end{split}\]
!pip install gurobipy
from gurobipy import * # This command imports the Gurobi functions and classes.
# create new model
m = Model('Western Clothing Company')
x1 = m.addVar(lb = 0, vtype = GRB.CONTINUOUS, name='x1')
x2 = m.addVar(lb = 0, vtype = GRB.CONTINUOUS, name='x2')
m.setObjective(((1500-x1)/24.6-12)*x1+((2700-x2)/63.8-9)*x2, GRB.MAXIMIZE)
m.addConstr(2*x1+2.7*x2<=6000, 'const1')
m.addConstr(3.6*x1+2.9*x2<=8500, 'const2')
m.addConstr(7.2*x1+8.5*x2<=15000, 'const3')
m.optimize()
#display optimal production plan
for v in m.getVars():
print(v.varName, v.x)
print('optimal total revenue:', m.objVal)
Requirement already satisfied: gurobipy in /usr/local/lib/python3.7/dist-packages (9.5.0)
Gurobi Optimizer version 9.5.0 build v9.5.0rc5 (linux64)
Thread count: 1 physical cores, 2 logical processors, using up to 2 threads
Optimize a model with 3 rows, 2 columns and 6 nonzeros
Model fingerprint: 0x390c5c3a
Model has 2 quadratic objective terms
Coefficient statistics:
Matrix range [2e+00, 8e+00]
Objective range [3e+01, 5e+01]
QObjective range [3e-02, 8e-02]
Bounds range [0e+00, 0e+00]
RHS range [6e+03, 2e+04]
Presolve time: 0.06s
Presolved: 3 rows, 2 columns, 6 nonzeros
Presolved model has 2 quadratic objective terms
Ordering time: 0.00s
Barrier statistics:
AA' NZ : 3.000e+00
Factor NZ : 6.000e+00
Factor Ops : 1.400e+01 (less than 1 second per iteration)
Threads : 1
Objective Residual
Iter Primal Dual Primal Dual Compl Time
0 -3.51759547e+05 6.37499504e+05 9.27e+02 0.00e+00 1.12e+06 0s
1 2.37322400e+04 6.43047090e+05 0.00e+00 0.00e+00 1.24e+05 0s
2 2.51606709e+04 6.37933000e+04 0.00e+00 0.00e+00 7.73e+03 0s
3 3.20504323e+04 4.25647642e+04 0.00e+00 0.00e+00 2.10e+03 0s
4 3.24417089e+04 3.27645350e+04 0.00e+00 0.00e+00 6.46e+01 0s
5 3.24592168e+04 3.24657001e+04 0.00e+00 1.42e-14 1.30e+00 0s
6 3.24592344e+04 3.24592408e+04 0.00e+00 0.00e+00 1.29e-03 0s
7 3.24592344e+04 3.24592344e+04 0.00e+00 2.84e-14 1.29e-06 0s
8 3.24592344e+04 3.24592344e+04 0.00e+00 0.00e+00 1.29e-09 0s
Barrier solved model in 8 iterations and 0.23 seconds (0.00 work units)
Optimal objective 3.24592344e+04
x1 602.3999999996935
x2 1062.8999999990601
optimal total revenue: 32459.234379539725
Problem Example: C10Q14¶
Solve the following nonlinear programming model:
\[\begin{split}
\begin{align}
&\text{max}\\
&\qquad z=\left(15000-\frac{9000}{x_1} \right)\\
&\qquad ~~+\left(24000-\frac{15000}{x_2} \right)\\
&\qquad ~~+\left(8100-\frac{5300}{x_3} \right)\\
&\qquad ~~+\left(12000-\frac{7600}{x_4} \right)\\
&\qquad ~~+\left(21000-\frac{12500}{x_4} \right)\\
&\text{s.t.}\\
&\qquad x_1+x_2+x_3+x_4+x_5 \le 15\\
&\qquad 355x_1+540x_2+290x_3+275x_4+490x_5 \le 6500\\
&\qquad x_i \ge 0,~\text{and integer}
\end{align}
\end{split}\]
from gurobipy import *
import numpy as np
m = Model('C10Q14')
m.params.NonConvex = 2
x = m.addVars(range(5), name='x', lb=1, vtype = GRB.INTEGER)
u = m.addVars(range(5), name='u', lb=0, vtype = GRB.CONTINUOUS)
m.setObjective((15000 -9000*u[0])+(24000-15000*u[1])+
(8100 -5300*u[2])+(12000 -7600*u[3])+
(21000-12500*u[4]), GRB.MAXIMIZE)
m.addConstr(x[0]+x[1]+x[2]+x[3]+x[4]<=15, 'const1')
m.addConstr(355*x[0]+540*x[1]+290*x[2]+275*x[3]+490*x[4]<=6500, 'const2')
m.addConstrs((x[i]*u[i] == 1 for i in range(5)), name='bilinear')
m.optimize()
for v in m.getVars():
print(v.varName, v.x)
print('optimal total revenue:', m.objVal)
Set parameter NonConvex to value 2
Gurobi Optimizer version 9.5.0 build v9.5.0rc5 (linux64)
Thread count: 1 physical cores, 2 logical processors, using up to 2 threads
Optimize a model with 2 rows, 10 columns and 10 nonzeros
Model fingerprint: 0xb4d06825
Model has 5 quadratic constraints
Variable types: 5 continuous, 5 integer (0 binary)
Coefficient statistics:
Matrix range [1e+00, 5e+02]
QMatrix range [1e+00, 1e+00]
Objective range [5e+03, 2e+04]
Bounds range [1e+00, 1e+00]
RHS range [2e+01, 6e+03]
QRHS range [1e+00, 1e+00]
Presolve time: 0.00s
Presolved: 22 rows, 10 columns, 50 nonzeros
Presolved model has 5 bilinear constraint(s)
Variable types: 5 continuous, 5 integer (0 binary)
Root relaxation: objective 7.218908e+04, 7 iterations, 0.00 seconds (0.00 work units)
Nodes | Current Node | Objective Bounds | Work
Expl Unexpl | Obj Depth IntInf | Incumbent BestBd Gap | It/Node Time
0 0 72189.0780 0 10 - 72189.0780 - - 0s
H 0 0 60400.000000 72189.0780 19.5% - 0s
0 0 71744.9856 0 10 60400.0000 71744.9856 18.8% - 0s
0 0 71058.7908 0 5 60400.0000 71058.7908 17.6% - 0s
H 0 0 63633.333333 71058.7908 11.7% - 0s
0 0 69577.4266 0 10 63633.3333 69577.4266 9.34% - 0s
0 0 69393.2479 0 5 63633.3333 69393.2479 9.05% - 0s
0 0 69352.0957 0 5 63633.3333 69352.0957 8.99% - 0s
H 0 0 63791.666667 69352.0957 8.72% - 0s
0 2 69352.0957 0 5 63791.6667 69352.0957 8.72% - 0s
H 22 10 64000.000000 67922.9290 6.13% 2.4 0s
Cutting planes:
MIR: 6
RLT: 6
Explored 70 nodes (193 simplex iterations) in 0.33 seconds (0.00 work units)
Thread count was 2 (of 2 available processors)
Solution count 4: 64000 63791.7 63633.3 60400
Optimal solution found (tolerance 1.00e-04)
Best objective 6.400000000000e+04, best bound 6.400000000000e+04, gap 0.0000%
x[0] 3.0
x[1] 4.0
x[2] 2.0
x[3] 3.0
x[4] 3.0
u[0] 0.3333333333333333
u[1] 0.25
u[2] 0.5
u[3] 0.3333333333333333
u[4] 0.3333333333333333
optimal total revenue: 64000.0
Facility Location Problem Example: C10Q19¶
Solve the following nonlinear programming model:
\[\begin{split}
\begin{align}
&\text{min}\\
&\qquad z=7000\sqrt{(1000-x)^2+(1250-y)^2}\\
&\qquad ~~+9000\sqrt{(1500-x)^2+(2700-y)^2}\\
&\qquad ~+11500\sqrt{(2000-x)^2+(700-y)^2}\\
&\qquad ~~+4300\sqrt{(2200-x)^2+(2000-y)^2}\\
\\
&\text{s.t.}\\
&\qquad x, y \ge 0
\end{align}
\end{split}\]
import gurobipy as gp
from gurobipy import GRB
# create new model
m = gp.Model('C10Q19')
m.params.NonConvex = 2
x = m.addVars(range(2), name='x', lb=0, vtype = GRB.CONTINUOUS)
u = m.addVars(range(4), name='u', vtype = GRB.CONTINUOUS)
v = m.addVars(range(4), name='v', vtype = GRB.CONTINUOUS)
m.setObjective(7000*v[0]+9000*v[1]+11500*v[2]+4300*v[3], GRB.MINIMIZE)
m.addConstr( (x[0]-1000)*(x[0]-1000)+(x[1]-1250)*(x[1]-1250) == u[0])
m.addConstr( (x[0]-1500)*(x[0]-1500)+(x[1]-2700)*(x[1]-2700) == u[1])
m.addConstr( (x[0]-2000)*(x[0]-2000)+(x[1]-700) *(x[1]-700) == u[2])
m.addConstr( (x[0]-2200)*(x[0]-2200)+(x[1]-2000)*(x[1]-2000) == u[3])
m.addConstrs((u[i] == v[i]*v[i] for i in range(4)), name='v_squared')
m.optimize()
#display optimal production plan
for w in m.getVars():
print(w.varName, w.x)
print('optimal total revenue:', m.objVal)
Set parameter NonConvex to value 2
Gurobi Optimizer version 9.5.0 build v9.5.0rc5 (linux64)
Thread count: 1 physical cores, 2 logical processors, using up to 2 threads
Optimize a model with 0 rows, 10 columns and 0 nonzeros
Model fingerprint: 0x087228d7
Model has 8 quadratic constraints
Coefficient statistics:
Matrix range [0e+00, 0e+00]
QMatrix range [1e+00, 1e+00]
QLMatrix range [1e+00, 5e+03]
Objective range [4e+03, 1e+04]
Bounds range [0e+00, 0e+00]
RHS range [0e+00, 0e+00]
QRHS range [3e+06, 1e+07]
Continuous model is non-convex -- solving as a MIP
Presolve time: 0.00s
Presolved: 22 rows, 13 columns, 50 nonzeros
Presolved model has 6 bilinear constraint(s)
Variable types: 13 continuous, 0 integer (0 binary)
Root relaxation: objective 0.000000e+00, 5 iterations, 0.00 seconds (0.00 work units)
Nodes | Current Node | Objective Bounds | Work
Expl Unexpl | Obj Depth IntInf | Incumbent BestBd Gap | It/Node Time
0 0 0.00000 0 4 - 0.00000 - - 0s
0 0 0.00000 0 4 - 0.00000 - - 0s
H 0 0 3.188442e+07 0.00000 100% - 0s
0 2 2955947.12 0 3 3.1884e+07 2955947.12 90.7% - 0s
H 38 28 2.904697e+07 2955947.12 89.8% 1.6 0s
H 99 51 2.902868e+07 2955947.12 89.8% 1.2 0s
H 194 56 2.901548e+07 8429605.93 70.9% 1.3 0s
* 781 107 54 2.896174e+07 2.4402e+07 15.7% 1.0 0s
* 1880 137 67 2.895822e+07 2.5464e+07 12.1% 1.0 0s
* 2608 121 66 2.894775e+07 2.8321e+07 2.17% 1.0 0s
* 2742 110 54 2.894588e+07 2.8500e+07 1.54% 1.0 0s
* 2866 113 62 2.894464e+07 2.8670e+07 0.95% 1.0 0s
* 3078 107 66 2.894464e+07 2.8858e+07 0.30% 1.0 0s
* 3079 105 66 2.894464e+07 2.8858e+07 0.30% 1.0 0s
Explored 3536 nodes (3724 simplex iterations) in 0.56 seconds (0.04 work units)
Thread count was 2 (of 2 available processors)
Solution count 10: 2.89446e+07 2.89446e+07 2.89446e+07 ... 2.9047e+07
Optimal solution found (tolerance 1.00e-04)
Best objective 2.894463512888e+07, best bound 2.894174284308e+07, gap 0.0100%
x[0] 1656.4309778346258
x[1] 1415.4868451162663
u[0] 458287.5245674583
u[1] 1674444.6958956604
u[2] 629961.0985260992
u[3] 637122.9100899571
v[0] 676.9693675251503
v[1] 1294.003360078937
v[2] 793.7008873156619
v[3] 798.1997933411493
optimal total revenue: 28944635.12888354
import numpy as np
x=1658.8;y=1416.7;
d=7000*np.sqrt((x-1000)**2+(y-1250)**2)+\
9000*np.sqrt((x-1500)**2+(y-2700)**2)+\
11500*np.sqrt((x-2000)**2+(y-700)**2)+\
4300*np.sqrt((x-2200)**2+(y-2000)**2)
print(x,y,d)
x=1665.4;y=1562.9;
d=7000*np.sqrt((x-1000)**2+(y-1250)**2)+\
9000*np.sqrt((x-1500)**2+(y-2700)**2)+\
11500*np.sqrt((x-2000)**2+(y-700)**2)+\
4300*np.sqrt((x-2200)**2+(y-2000)**2)
print(x,y,d)
1658.8 1416.7 28944633.637902677
1665.4 1562.9 29101303.298960045