Coverage for /usr/lib/python3/dist-packages/scipy/optimize/_linprog.py: 16%

103 statements  

« prev     ^ index     » next       coverage.py v7.9.1, created at 2025-06-14 15:55 +0200

1""" 

2A top-level linear programming interface. 

3 

4.. versionadded:: 0.15.0 

5 

6Functions 

7--------- 

8.. autosummary:: 

9 :toctree: generated/ 

10 

11 linprog 

12 linprog_verbose_callback 

13 linprog_terse_callback 

14 

15""" 

16 

17import numpy as np 

18 

19from ._optimize import OptimizeResult, OptimizeWarning 

20from warnings import warn 

21from ._linprog_highs import _linprog_highs 

22from ._linprog_ip import _linprog_ip 

23from ._linprog_simplex import _linprog_simplex 

24from ._linprog_rs import _linprog_rs 

25from ._linprog_doc import (_linprog_highs_doc, _linprog_ip_doc, 

26 _linprog_rs_doc, _linprog_simplex_doc, 

27 _linprog_highs_ipm_doc, _linprog_highs_ds_doc) 

28from ._linprog_util import ( 

29 _parse_linprog, _presolve, _get_Abc, _LPProblem, _autoscale, 

30 _postsolve, _check_result, _display_summary) 

31from copy import deepcopy 

32 

33__all__ = ['linprog', 'linprog_verbose_callback', 'linprog_terse_callback'] 

34 

35__docformat__ = "restructuredtext en" 

36 

37LINPROG_METHODS = ['simplex', 'revised simplex', 'interior-point', 'highs', 'highs-ds', 'highs-ipm'] 

38 

39 

40def linprog_verbose_callback(res): 

41 """ 

42 A sample callback function demonstrating the linprog callback interface. 

43 This callback produces detailed output to sys.stdout before each iteration 

44 and after the final iteration of the simplex algorithm. 

45 

46 Parameters 

47 ---------- 

48 res : A `scipy.optimize.OptimizeResult` consisting of the following fields: 

49 

50 x : 1-D array 

51 The independent variable vector which optimizes the linear 

52 programming problem. 

53 fun : float 

54 Value of the objective function. 

55 success : bool 

56 True if the algorithm succeeded in finding an optimal solution. 

57 slack : 1-D array 

58 The values of the slack variables. Each slack variable corresponds 

59 to an inequality constraint. If the slack is zero, then the 

60 corresponding constraint is active. 

61 con : 1-D array 

62 The (nominally zero) residuals of the equality constraints, that is, 

63 ``b - A_eq @ x`` 

64 phase : int 

65 The phase of the optimization being executed. In phase 1 a basic 

66 feasible solution is sought and the T has an additional row 

67 representing an alternate objective function. 

68 status : int 

69 An integer representing the exit status of the optimization:: 

70 

71 0 : Optimization terminated successfully 

72 1 : Iteration limit reached 

73 2 : Problem appears to be infeasible 

74 3 : Problem appears to be unbounded 

75 4 : Serious numerical difficulties encountered 

76 

77 nit : int 

78 The number of iterations performed. 

79 message : str 

80 A string descriptor of the exit status of the optimization. 

81 """ 

82 x = res['x'] 

83 fun = res['fun'] 

84 phase = res['phase'] 

85 status = res['status'] 

86 nit = res['nit'] 

87 message = res['message'] 

88 complete = res['complete'] 

89 

90 saved_printoptions = np.get_printoptions() 

91 np.set_printoptions(linewidth=500, 

92 formatter={'float': lambda x: f"{x: 12.4f}"}) 

93 if status: 

94 print('--------- Simplex Early Exit -------\n') 

95 print(f'The simplex method exited early with status {status:d}') 

96 print(message) 

97 elif complete: 

98 print('--------- Simplex Complete --------\n') 

99 print(f'Iterations required: {nit}') 

100 else: 

101 print(f'--------- Iteration {nit:d} ---------\n') 

102 

103 if nit > 0: 

104 if phase == 1: 

105 print('Current Pseudo-Objective Value:') 

106 else: 

107 print('Current Objective Value:') 

108 print('f = ', fun) 

109 print() 

110 print('Current Solution Vector:') 

111 print('x = ', x) 

112 print() 

113 

114 np.set_printoptions(**saved_printoptions) 

115 

116 

117def linprog_terse_callback(res): 

118 """ 

119 A sample callback function demonstrating the linprog callback interface. 

120 This callback produces brief output to sys.stdout before each iteration 

121 and after the final iteration of the simplex algorithm. 

122 

123 Parameters 

124 ---------- 

125 res : A `scipy.optimize.OptimizeResult` consisting of the following fields: 

126 

127 x : 1-D array 

128 The independent variable vector which optimizes the linear 

129 programming problem. 

130 fun : float 

131 Value of the objective function. 

132 success : bool 

133 True if the algorithm succeeded in finding an optimal solution. 

134 slack : 1-D array 

135 The values of the slack variables. Each slack variable corresponds 

136 to an inequality constraint. If the slack is zero, then the 

137 corresponding constraint is active. 

138 con : 1-D array 

139 The (nominally zero) residuals of the equality constraints, that is, 

140 ``b - A_eq @ x``. 

141 phase : int 

142 The phase of the optimization being executed. In phase 1 a basic 

143 feasible solution is sought and the T has an additional row 

144 representing an alternate objective function. 

145 status : int 

146 An integer representing the exit status of the optimization:: 

147 

148 0 : Optimization terminated successfully 

149 1 : Iteration limit reached 

150 2 : Problem appears to be infeasible 

151 3 : Problem appears to be unbounded 

152 4 : Serious numerical difficulties encountered 

153 

154 nit : int 

155 The number of iterations performed. 

156 message : str 

157 A string descriptor of the exit status of the optimization. 

158 """ 

159 nit = res['nit'] 

160 x = res['x'] 

161 

162 if nit == 0: 

163 print("Iter: X:") 

164 print(f"{nit: <5d} ", end="") 

165 print(x) 

166 

167 

168def linprog(c, A_ub=None, b_ub=None, A_eq=None, b_eq=None, 

169 bounds=None, method='highs', callback=None, 

170 options=None, x0=None, integrality=None): 

171 r""" 

172 Linear programming: minimize a linear objective function subject to linear 

173 equality and inequality constraints. 

174 

175 Linear programming solves problems of the following form: 

176 

177 .. math:: 

178 

179 \min_x \ & c^T x \\ 

180 \mbox{such that} \ & A_{ub} x \leq b_{ub},\\ 

181 & A_{eq} x = b_{eq},\\ 

182 & l \leq x \leq u , 

183 

184 where :math:`x` is a vector of decision variables; :math:`c`, 

185 :math:`b_{ub}`, :math:`b_{eq}`, :math:`l`, and :math:`u` are vectors; and 

186 :math:`A_{ub}` and :math:`A_{eq}` are matrices. 

187 

188 Alternatively, that's: 

189 

190 - minimize :: 

191 

192 c @ x 

193 

194 - such that :: 

195 

196 A_ub @ x <= b_ub 

197 A_eq @ x == b_eq 

198 lb <= x <= ub 

199 

200 Note that by default ``lb = 0`` and ``ub = None``. Other bounds can be 

201 specified with ``bounds``. 

202 

203 Parameters 

204 ---------- 

205 c : 1-D array 

206 The coefficients of the linear objective function to be minimized. 

207 A_ub : 2-D array, optional 

208 The inequality constraint matrix. Each row of ``A_ub`` specifies the 

209 coefficients of a linear inequality constraint on ``x``. 

210 b_ub : 1-D array, optional 

211 The inequality constraint vector. Each element represents an 

212 upper bound on the corresponding value of ``A_ub @ x``. 

213 A_eq : 2-D array, optional 

214 The equality constraint matrix. Each row of ``A_eq`` specifies the 

215 coefficients of a linear equality constraint on ``x``. 

216 b_eq : 1-D array, optional 

217 The equality constraint vector. Each element of ``A_eq @ x`` must equal 

218 the corresponding element of ``b_eq``. 

219 bounds : sequence, optional 

220 A sequence of ``(min, max)`` pairs for each element in ``x``, defining 

221 the minimum and maximum values of that decision variable. 

222 If a single tuple ``(min, max)`` is provided, then ``min`` and ``max`` 

223 will serve as bounds for all decision variables. 

224 Use ``None`` to indicate that there is no bound. For instance, the 

225 default bound ``(0, None)`` means that all decision variables are 

226 non-negative, and the pair ``(None, None)`` means no bounds at all, 

227 i.e. all variables are allowed to be any real. 

228 method : str, optional 

229 The algorithm used to solve the standard form problem. 

230 :ref:`'highs' <optimize.linprog-highs>` (default), 

231 :ref:`'highs-ds' <optimize.linprog-highs-ds>`, 

232 :ref:`'highs-ipm' <optimize.linprog-highs-ipm>`, 

233 :ref:`'interior-point' <optimize.linprog-interior-point>` (legacy), 

234 :ref:`'revised simplex' <optimize.linprog-revised_simplex>` (legacy), 

235 and 

236 :ref:`'simplex' <optimize.linprog-simplex>` (legacy) are supported. 

237 The legacy methods are deprecated and will be removed in SciPy 1.11.0. 

238 callback : callable, optional 

239 If a callback function is provided, it will be called at least once per 

240 iteration of the algorithm. The callback function must accept a single 

241 `scipy.optimize.OptimizeResult` consisting of the following fields: 

242 

243 x : 1-D array 

244 The current solution vector. 

245 fun : float 

246 The current value of the objective function ``c @ x``. 

247 success : bool 

248 ``True`` when the algorithm has completed successfully. 

249 slack : 1-D array 

250 The (nominally positive) values of the slack, 

251 ``b_ub - A_ub @ x``. 

252 con : 1-D array 

253 The (nominally zero) residuals of the equality constraints, 

254 ``b_eq - A_eq @ x``. 

255 phase : int 

256 The phase of the algorithm being executed. 

257 status : int 

258 An integer representing the status of the algorithm. 

259 

260 ``0`` : Optimization proceeding nominally. 

261 

262 ``1`` : Iteration limit reached. 

263 

264 ``2`` : Problem appears to be infeasible. 

265 

266 ``3`` : Problem appears to be unbounded. 

267 

268 ``4`` : Numerical difficulties encountered. 

269 

270 nit : int 

271 The current iteration number. 

272 message : str 

273 A string descriptor of the algorithm status. 

274 

275 Callback functions are not currently supported by the HiGHS methods. 

276 

277 options : dict, optional 

278 A dictionary of solver options. All methods accept the following 

279 options: 

280 

281 maxiter : int 

282 Maximum number of iterations to perform. 

283 Default: see method-specific documentation. 

284 disp : bool 

285 Set to ``True`` to print convergence messages. 

286 Default: ``False``. 

287 presolve : bool 

288 Set to ``False`` to disable automatic presolve. 

289 Default: ``True``. 

290 

291 All methods except the HiGHS solvers also accept: 

292 

293 tol : float 

294 A tolerance which determines when a residual is "close enough" to 

295 zero to be considered exactly zero. 

296 autoscale : bool 

297 Set to ``True`` to automatically perform equilibration. 

298 Consider using this option if the numerical values in the 

299 constraints are separated by several orders of magnitude. 

300 Default: ``False``. 

301 rr : bool 

302 Set to ``False`` to disable automatic redundancy removal. 

303 Default: ``True``. 

304 rr_method : string 

305 Method used to identify and remove redundant rows from the 

306 equality constraint matrix after presolve. For problems with 

307 dense input, the available methods for redundancy removal are: 

308 

309 "SVD": 

310 Repeatedly performs singular value decomposition on 

311 the matrix, detecting redundant rows based on nonzeros 

312 in the left singular vectors that correspond with 

313 zero singular values. May be fast when the matrix is 

314 nearly full rank. 

315 "pivot": 

316 Uses the algorithm presented in [5]_ to identify 

317 redundant rows. 

318 "ID": 

319 Uses a randomized interpolative decomposition. 

320 Identifies columns of the matrix transpose not used in 

321 a full-rank interpolative decomposition of the matrix. 

322 None: 

323 Uses "svd" if the matrix is nearly full rank, that is, 

324 the difference between the matrix rank and the number 

325 of rows is less than five. If not, uses "pivot". The 

326 behavior of this default is subject to change without 

327 prior notice. 

328 

329 Default: None. 

330 For problems with sparse input, this option is ignored, and the 

331 pivot-based algorithm presented in [5]_ is used. 

332 

333 For method-specific options, see 

334 :func:`show_options('linprog') <show_options>`. 

335 

336 x0 : 1-D array, optional 

337 Guess values of the decision variables, which will be refined by 

338 the optimization algorithm. This argument is currently used only by the 

339 'revised simplex' method, and can only be used if `x0` represents a 

340 basic feasible solution. 

341 

342 integrality : 1-D array or int, optional 

343 Indicates the type of integrality constraint on each decision variable. 

344 

345 ``0`` : Continuous variable; no integrality constraint. 

346 

347 ``1`` : Integer variable; decision variable must be an integer 

348 within `bounds`. 

349 

350 ``2`` : Semi-continuous variable; decision variable must be within 

351 `bounds` or take value ``0``. 

352 

353 ``3`` : Semi-integer variable; decision variable must be an integer 

354 within `bounds` or take value ``0``. 

355 

356 By default, all variables are continuous. 

357 

358 For mixed integrality constraints, supply an array of shape `c.shape`. 

359 To infer a constraint on each decision variable from shorter inputs, 

360 the argument will be broadcasted to `c.shape` using `np.broadcast_to`. 

361 

362 This argument is currently used only by the ``'highs'`` method and 

363 ignored otherwise. 

364 

365 Returns 

366 ------- 

367 res : OptimizeResult 

368 A :class:`scipy.optimize.OptimizeResult` consisting of the fields 

369 below. Note that the return types of the fields may depend on whether 

370 the optimization was successful, therefore it is recommended to check 

371 `OptimizeResult.status` before relying on the other fields: 

372 

373 x : 1-D array 

374 The values of the decision variables that minimizes the 

375 objective function while satisfying the constraints. 

376 fun : float 

377 The optimal value of the objective function ``c @ x``. 

378 slack : 1-D array 

379 The (nominally positive) values of the slack variables, 

380 ``b_ub - A_ub @ x``. 

381 con : 1-D array 

382 The (nominally zero) residuals of the equality constraints, 

383 ``b_eq - A_eq @ x``. 

384 success : bool 

385 ``True`` when the algorithm succeeds in finding an optimal 

386 solution. 

387 status : int 

388 An integer representing the exit status of the algorithm. 

389 

390 ``0`` : Optimization terminated successfully. 

391 

392 ``1`` : Iteration limit reached. 

393 

394 ``2`` : Problem appears to be infeasible. 

395 

396 ``3`` : Problem appears to be unbounded. 

397 

398 ``4`` : Numerical difficulties encountered. 

399 

400 nit : int 

401 The total number of iterations performed in all phases. 

402 message : str 

403 A string descriptor of the exit status of the algorithm. 

404 

405 See Also 

406 -------- 

407 show_options : Additional options accepted by the solvers. 

408 

409 Notes 

410 ----- 

411 This section describes the available solvers that can be selected by the 

412 'method' parameter. 

413 

414 `'highs-ds'` and 

415 `'highs-ipm'` are interfaces to the 

416 HiGHS simplex and interior-point method solvers [13]_, respectively. 

417 `'highs'` (default) chooses between 

418 the two automatically. These are the fastest linear 

419 programming solvers in SciPy, especially for large, sparse problems; 

420 which of these two is faster is problem-dependent. 

421 The other solvers (`'interior-point'`, `'revised simplex'`, and 

422 `'simplex'`) are legacy methods and will be removed in SciPy 1.11.0. 

423 

424 Method *highs-ds* is a wrapper of the C++ high performance dual 

425 revised simplex implementation (HSOL) [13]_, [14]_. Method *highs-ipm* 

426 is a wrapper of a C++ implementation of an **i**\ nterior-\ **p**\ oint 

427 **m**\ ethod [13]_; it features a crossover routine, so it is as accurate 

428 as a simplex solver. Method *highs* chooses between the two automatically. 

429 For new code involving `linprog`, we recommend explicitly choosing one of 

430 these three method values. 

431 

432 .. versionadded:: 1.6.0 

433 

434 Method *interior-point* uses the primal-dual path following algorithm 

435 as outlined in [4]_. This algorithm supports sparse constraint matrices and 

436 is typically faster than the simplex methods, especially for large, sparse 

437 problems. Note, however, that the solution returned may be slightly less 

438 accurate than those of the simplex methods and will not, in general, 

439 correspond with a vertex of the polytope defined by the constraints. 

440 

441 .. versionadded:: 1.0.0 

442 

443 Method *revised simplex* uses the revised simplex method as described in 

444 [9]_, except that a factorization [11]_ of the basis matrix, rather than 

445 its inverse, is efficiently maintained and used to solve the linear systems 

446 at each iteration of the algorithm. 

447 

448 .. versionadded:: 1.3.0 

449 

450 Method *simplex* uses a traditional, full-tableau implementation of 

451 Dantzig's simplex algorithm [1]_, [2]_ (*not* the 

452 Nelder-Mead simplex). This algorithm is included for backwards 

453 compatibility and educational purposes. 

454 

455 .. versionadded:: 0.15.0 

456 

457 Before applying *interior-point*, *revised simplex*, or *simplex*, 

458 a presolve procedure based on [8]_ attempts 

459 to identify trivial infeasibilities, trivial unboundedness, and potential 

460 problem simplifications. Specifically, it checks for: 

461 

462 - rows of zeros in ``A_eq`` or ``A_ub``, representing trivial constraints; 

463 - columns of zeros in ``A_eq`` `and` ``A_ub``, representing unconstrained 

464 variables; 

465 - column singletons in ``A_eq``, representing fixed variables; and 

466 - column singletons in ``A_ub``, representing simple bounds. 

467 

468 If presolve reveals that the problem is unbounded (e.g. an unconstrained 

469 and unbounded variable has negative cost) or infeasible (e.g., a row of 

470 zeros in ``A_eq`` corresponds with a nonzero in ``b_eq``), the solver 

471 terminates with the appropriate status code. Note that presolve terminates 

472 as soon as any sign of unboundedness is detected; consequently, a problem 

473 may be reported as unbounded when in reality the problem is infeasible 

474 (but infeasibility has not been detected yet). Therefore, if it is 

475 important to know whether the problem is actually infeasible, solve the 

476 problem again with option ``presolve=False``. 

477 

478 If neither infeasibility nor unboundedness are detected in a single pass 

479 of the presolve, bounds are tightened where possible and fixed 

480 variables are removed from the problem. Then, linearly dependent rows 

481 of the ``A_eq`` matrix are removed, (unless they represent an 

482 infeasibility) to avoid numerical difficulties in the primary solve 

483 routine. Note that rows that are nearly linearly dependent (within a 

484 prescribed tolerance) may also be removed, which can change the optimal 

485 solution in rare cases. If this is a concern, eliminate redundancy from 

486 your problem formulation and run with option ``rr=False`` or 

487 ``presolve=False``. 

488 

489 Several potential improvements can be made here: additional presolve 

490 checks outlined in [8]_ should be implemented, the presolve routine should 

491 be run multiple times (until no further simplifications can be made), and 

492 more of the efficiency improvements from [5]_ should be implemented in the 

493 redundancy removal routines. 

494 

495 After presolve, the problem is transformed to standard form by converting 

496 the (tightened) simple bounds to upper bound constraints, introducing 

497 non-negative slack variables for inequality constraints, and expressing 

498 unbounded variables as the difference between two non-negative variables. 

499 Optionally, the problem is automatically scaled via equilibration [12]_. 

500 The selected algorithm solves the standard form problem, and a 

501 postprocessing routine converts the result to a solution to the original 

502 problem. 

503 

504 References 

505 ---------- 

506 .. [1] Dantzig, George B., Linear programming and extensions. Rand 

507 Corporation Research Study Princeton Univ. Press, Princeton, NJ, 

508 1963 

509 .. [2] Hillier, S.H. and Lieberman, G.J. (1995), "Introduction to 

510 Mathematical Programming", McGraw-Hill, Chapter 4. 

511 .. [3] Bland, Robert G. New finite pivoting rules for the simplex method. 

512 Mathematics of Operations Research (2), 1977: pp. 103-107. 

513 .. [4] Andersen, Erling D., and Knud D. Andersen. "The MOSEK interior point 

514 optimizer for linear programming: an implementation of the 

515 homogeneous algorithm." High performance optimization. Springer US, 

516 2000. 197-232. 

517 .. [5] Andersen, Erling D. "Finding all linearly dependent rows in 

518 large-scale linear programming." Optimization Methods and Software 

519 6.3 (1995): 219-227. 

520 .. [6] Freund, Robert M. "Primal-Dual Interior-Point Methods for Linear 

521 Programming based on Newton's Method." Unpublished Course Notes, 

522 March 2004. Available 2/25/2017 at 

523 https://ocw.mit.edu/courses/sloan-school-of-management/15-084j-nonlinear-programming-spring-2004/lecture-notes/lec14_int_pt_mthd.pdf 

524 .. [7] Fourer, Robert. "Solving Linear Programs by Interior-Point Methods." 

525 Unpublished Course Notes, August 26, 2005. Available 2/25/2017 at 

526 http://www.4er.org/CourseNotes/Book%20B/B-III.pdf 

527 .. [8] Andersen, Erling D., and Knud D. Andersen. "Presolving in linear 

528 programming." Mathematical Programming 71.2 (1995): 221-245. 

529 .. [9] Bertsimas, Dimitris, and J. Tsitsiklis. "Introduction to linear 

530 programming." Athena Scientific 1 (1997): 997. 

531 .. [10] Andersen, Erling D., et al. Implementation of interior point 

532 methods for large scale linear programming. HEC/Universite de 

533 Geneve, 1996. 

534 .. [11] Bartels, Richard H. "A stabilization of the simplex method." 

535 Journal in Numerische Mathematik 16.5 (1971): 414-434. 

536 .. [12] Tomlin, J. A. "On scaling linear programming problems." 

537 Mathematical Programming Study 4 (1975): 146-166. 

538 .. [13] Huangfu, Q., Galabova, I., Feldmeier, M., and Hall, J. A. J. 

539 "HiGHS - high performance software for linear optimization." 

540 https://highs.dev/ 

541 .. [14] Huangfu, Q. and Hall, J. A. J. "Parallelizing the dual revised 

542 simplex method." Mathematical Programming Computation, 10 (1), 

543 119-142, 2018. DOI: 10.1007/s12532-017-0130-5 

544 

545 Examples 

546 -------- 

547 Consider the following problem: 

548 

549 .. math:: 

550 

551 \min_{x_0, x_1} \ -x_0 + 4x_1 & \\ 

552 \mbox{such that} \ -3x_0 + x_1 & \leq 6,\\ 

553 -x_0 - 2x_1 & \geq -4,\\ 

554 x_1 & \geq -3. 

555 

556 The problem is not presented in the form accepted by `linprog`. This is 

557 easily remedied by converting the "greater than" inequality 

558 constraint to a "less than" inequality constraint by 

559 multiplying both sides by a factor of :math:`-1`. Note also that the last 

560 constraint is really the simple bound :math:`-3 \leq x_1 \leq \infty`. 

561 Finally, since there are no bounds on :math:`x_0`, we must explicitly 

562 specify the bounds :math:`-\infty \leq x_0 \leq \infty`, as the 

563 default is for variables to be non-negative. After collecting coeffecients 

564 into arrays and tuples, the input for this problem is: 

565 

566 >>> from scipy.optimize import linprog 

567 >>> c = [-1, 4] 

568 >>> A = [[-3, 1], [1, 2]] 

569 >>> b = [6, 4] 

570 >>> x0_bounds = (None, None) 

571 >>> x1_bounds = (-3, None) 

572 >>> res = linprog(c, A_ub=A, b_ub=b, bounds=[x0_bounds, x1_bounds]) 

573 >>> res.fun 

574 -22.0 

575 >>> res.x 

576 array([10., -3.]) 

577 >>> res.message 

578 'Optimization terminated successfully. (HiGHS Status 7: Optimal)' 

579 

580 The marginals (AKA dual values / shadow prices / Lagrange multipliers) 

581 and residuals (slacks) are also available. 

582 

583 >>> res.ineqlin 

584 residual: [ 3.900e+01 0.000e+00] 

585 marginals: [-0.000e+00 -1.000e+00] 

586 

587 For example, because the marginal associated with the second inequality 

588 constraint is -1, we expect the optimal value of the objective function 

589 to decrease by ``eps`` if we add a small amount ``eps`` to the right hand 

590 side of the second inequality constraint: 

591 

592 >>> eps = 0.05 

593 >>> b[1] += eps 

594 >>> linprog(c, A_ub=A, b_ub=b, bounds=[x0_bounds, x1_bounds]).fun 

595 -22.05 

596 

597 Also, because the residual on the first inequality constraint is 39, we 

598 can decrease the right hand side of the first constraint by 39 without 

599 affecting the optimal solution. 

600 

601 >>> b = [6, 4] # reset to original values 

602 >>> b[0] -= 39 

603 >>> linprog(c, A_ub=A, b_ub=b, bounds=[x0_bounds, x1_bounds]).fun 

604 -22.0 

605 

606 """ 

607 

608 meth = method.lower() 

609 methods = {"highs", "highs-ds", "highs-ipm", 

610 "simplex", "revised simplex", "interior-point"} 

611 

612 if meth not in methods: 

613 raise ValueError(f"Unknown solver '{method}'") 

614 

615 if x0 is not None and meth != "revised simplex": 

616 warning_message = "x0 is used only when method is 'revised simplex'. " 

617 warn(warning_message, OptimizeWarning) 

618 

619 if np.any(integrality) and not meth == "highs": 

620 integrality = None 

621 warning_message = ("Only `method='highs'` supports integer " 

622 "constraints. Ignoring `integrality`.") 

623 warn(warning_message, OptimizeWarning) 

624 elif np.any(integrality): 

625 integrality = np.broadcast_to(integrality, np.shape(c)) 

626 

627 lp = _LPProblem(c, A_ub, b_ub, A_eq, b_eq, bounds, x0, integrality) 

628 lp, solver_options = _parse_linprog(lp, options, meth) 

629 tol = solver_options.get('tol', 1e-9) 

630 

631 # Give unmodified problem to HiGHS 

632 if meth.startswith('highs'): 

633 if callback is not None: 

634 raise NotImplementedError("HiGHS solvers do not support the " 

635 "callback interface.") 

636 highs_solvers = {'highs-ipm': 'ipm', 'highs-ds': 'simplex', 

637 'highs': None} 

638 

639 sol = _linprog_highs(lp, solver=highs_solvers[meth], 

640 **solver_options) 

641 sol['status'], sol['message'] = ( 

642 _check_result(sol['x'], sol['fun'], sol['status'], sol['slack'], 

643 sol['con'], lp.bounds, tol, sol['message'], 

644 integrality)) 

645 sol['success'] = sol['status'] == 0 

646 return OptimizeResult(sol) 

647 

648 warn(f"`method='{meth}'` is deprecated and will be removed in SciPy " 

649 "1.11.0. Please use one of the HiGHS solvers (e.g. " 

650 "`method='highs'`) in new code.", DeprecationWarning, stacklevel=2) 

651 

652 iteration = 0 

653 complete = False # will become True if solved in presolve 

654 undo = [] 

655 

656 # Keep the original arrays to calculate slack/residuals for original 

657 # problem. 

658 lp_o = deepcopy(lp) 

659 

660 # Solve trivial problem, eliminate variables, tighten bounds, etc. 

661 rr_method = solver_options.pop('rr_method', None) # need to pop these; 

662 rr = solver_options.pop('rr', True) # they're not passed to methods 

663 c0 = 0 # we might get a constant term in the objective 

664 if solver_options.pop('presolve', True): 

665 (lp, c0, x, undo, complete, status, message) = _presolve(lp, rr, 

666 rr_method, 

667 tol) 

668 

669 C, b_scale = 1, 1 # for trivial unscaling if autoscale is not used 

670 postsolve_args = (lp_o._replace(bounds=lp.bounds), undo, C, b_scale) 

671 

672 if not complete: 

673 A, b, c, c0, x0 = _get_Abc(lp, c0) 

674 if solver_options.pop('autoscale', False): 

675 A, b, c, x0, C, b_scale = _autoscale(A, b, c, x0) 

676 postsolve_args = postsolve_args[:-2] + (C, b_scale) 

677 

678 if meth == 'simplex': 

679 x, status, message, iteration = _linprog_simplex( 

680 c, c0=c0, A=A, b=b, callback=callback, 

681 postsolve_args=postsolve_args, **solver_options) 

682 elif meth == 'interior-point': 

683 x, status, message, iteration = _linprog_ip( 

684 c, c0=c0, A=A, b=b, callback=callback, 

685 postsolve_args=postsolve_args, **solver_options) 

686 elif meth == 'revised simplex': 

687 x, status, message, iteration = _linprog_rs( 

688 c, c0=c0, A=A, b=b, x0=x0, callback=callback, 

689 postsolve_args=postsolve_args, **solver_options) 

690 

691 # Eliminate artificial variables, re-introduce presolved variables, etc. 

692 disp = solver_options.get('disp', False) 

693 

694 x, fun, slack, con = _postsolve(x, postsolve_args, complete) 

695 

696 status, message = _check_result(x, fun, status, slack, con, lp_o.bounds, 

697 tol, message, integrality) 

698 

699 if disp: 

700 _display_summary(message, status, fun, iteration) 

701 

702 sol = { 

703 'x': x, 

704 'fun': fun, 

705 'slack': slack, 

706 'con': con, 

707 'status': status, 

708 'message': message, 

709 'nit': iteration, 

710 'success': status == 0} 

711 

712 return OptimizeResult(sol)