nAG数値計算ライブラリ
> 最適化アルゴリズムExample集
> 密な非線形計画問題の解法 (nlp2_solve)
非線形計画法-NLP(密)
このExampleは、非線形計画問題(NLP)のHockとSchittkowskiの問題71を、nAGライブラリのnlp2_solveを用いて解いています。
目的関数:
| タスク | 式 |
|---|---|
| minimize | \(x_0x_3(x_0 + x_1 + x_2) + x_2\) |
決定変数:
| 変数 | 範囲 |
|---|---|
| \(x_0\) | \(1 \leq x_0 \leq 5\) |
| \(x_1\) | \(1 \leq x_1 \leq 5\) |
| \(x_2\) | \(1 \leq x_2 \leq 5\) |
| \(x_3\) | \(1 \leq x_3 \leq 5\) |
制約条件:
| 制約 | 式 |
|---|---|
| 非線形制約1 | \(x_0^2 + x_1^2 + x_2^2 + x_3^2 \leq 40\) |
| 非線形制約2 | \(x_0x_1x_2x_3 \geq 25\) |
| 線形制約 | \(1 \leq x_0 + x_1 + x_2 + x_3 \leq 5\) |
コードでは、決定変数の初期値を \(x=[1, 5,
5, 1]\) と設定しています。また、変数の上下限制約については
bl と bu で定義されており、非線形制約の数は
ncnln=2 で指定されています。目的関数と非線形制約はそれぞれ
cb_objfun と cb_confun
で定義されています。
Exampleの実行コマンド:
python -m naginterfaces.library.examples.opt.nlp2_solve_ex
ソースコード表示コマンド:
python -c "import inspect; from naginterfaces.library.examples.opt import nlp2_solve_ex; print(''.join(inspect.getsourcelines(nlp2_solve_ex)[0]))"
出力結果例:
naginterfaces.library.opt.nlp2_solve Python Example Results.
Solve Hock and Schittkowski Problem 71.
Final objective value is 1.7014017e+01
マニュアル:
ソース:
#!/usr/bin/env python3
"``naginterfaces.library.opt.nlp2_solve`` Python Example."
# nAG Copyright 2017-2019.
# pylint: disable=invalid-name,too-many-locals
import numpy as np
from naginterfaces.library import opt
def main():
"""
Example for :func:`naginterfaces.library.opt.nlp2_solve`.
Dense NLP.
>>> main()
naginterfaces.library.opt.nlp2_solve Python Example Results.
Solve Hock and Schittkowski Problem 71.
Final objective value is 1.7014017e+01
"""
print(
'naginterfaces.library.opt.nlp2_solve Python Example Results.'
)
print('Solve Hock and Schittkowski Problem 71.')
def cb_confun(mode, needc, x, cjac, _nstate):
"""The nonlinear constraints."""
ccon = np.empty(max(1, len(needc)))
if needc[0] > 0:
if mode in [0, 2]:
ccon[0] = (x[0]**2 + x[1]**2 + x[2]**2 + x[3]**2)
if mode == 2:
cjac[0, :] = 2*x
if needc[1] > 0:
if mode in [0, 2]:
ccon[1] = x[0]*x[1]*x[2]*x[3]
if mode == 2:
cjac[1, :] = [
x[1]*x[2]*x[3],
x[0]*x[2]*x[3],
x[0]*x[1]*x[3],
x[0]*x[1]*x[2],
]
return 0, ccon, cjac
def cb_objfun(mode, x, grad, _nstate):
"""The objective function."""
if mode in [0, 2]:
objf = x[0]*x[3]*(x[0] + x[1] + x[2]) + x[2]
else:
objf = 0.
if mode == 2:
grad[:] = [
x[3]*(2*x[0] + x[1] + x[2]),
x[0]*x[3],
x[0]*x[3] + 1.0,
x[0]*(x[0] + x[1] + x[2]),
]
return 0, objf, grad
# Initialize the solver:
comm = opt.nlp2_init()
# The initial guess:
x = [1., 5., 5., 1.]
# The linear constraints:
a = np.array([[1.]*len(x)])
# There are two nonlinear constraints defined by cb_confun:
ncnln = 2
# The bounds:
bl = [1., 1., 1., 1., -1.0E+25, -1.0E+25, 25.]
bu = [5., 5., 5., 5., 20., 40., 1.0E+25]
n = len(x)
istate = np.zeros(n + 1 + ncnln, dtype=int)
ccon = np.zeros(max(1, ncnln))
cjac = np.zeros((ncnln, n))
clamda = np.zeros(n + 1 + ncnln)
h = np.zeros((n, n))
objf = opt.nlp2_solve(
a, bl, bu, cb_objfun,
istate, ccon, cjac, clamda, h, x, comm,
confun=cb_confun,
).objf
print('Final objective value is {:.7e}'.format(objf))
if __name__ == '__main__':
import doctest
import sys
sys.exit(
doctest.testmod(
None, verbose=True, report=False,
optionflags=doctest.REPORT_NDIFF,
).failed
)