nAG数値計算ライブラリ
> 最適化アルゴリズムExample集
> 疎な非線形最適化問題の逐次二次計画法による解法
逐次二次計画法-SQP（疎）
このExampleは、スパースな非線形計画問題(NLP)を逐次二次計画法(SQP)で解いています。2次の目的関数、線形制約、2つの非線形制約を含む問題を扱っています。目的関数は、ユーザコールバックの使用方法を示すために非線形関数としてコード化されています。
目的関数:
| タスク | 式 |
|---|---|
| minimize | \(f(x) = (x_0 + x_1 + x_2)^2 + 3x_2 + 5x_3 + \cos(0.01x_0) - 1\) |
決定変数:
| 変数 | 範囲 |
|---|---|
| \(x_0\) | \((-\infty, \infty)\) |
| \(x_1\) | \((-\infty, \infty)\) |
| \(x_2\) | \([0, \infty)\) |
| \(x_3\) | \([0, \infty)\) |
制約条件:
| 制約 | 式 |
|---|---|
| 線形制約 | \(2x_0 + 4x_1 \geq 0\) |
| 非線形制約1 | \(x_0^2 + x_1^2 + x_2 = 2\) |
| 非線形制約2 | \(x_1^4 + x_3 = 4\) |
目的関数は2次形式と非線形項の和で構成されています。\(x_0, x_1, x_2\) の2次形式に加えて、\(x_2, x_3\) の線形項、および \(x_0\) の非線形項 \(\cos(0.01x_0)\) が含まれています。
制約条件は1つの線形不等式制約と2つの非線形等式制約から成ります。線形制約は \(x_0, x_1\) に関する不等式で表現されています。非線形制約1は \(x_0, x_1, x_2\) に関する2次形式と定数項からなる等式です。非線形制約2は \(x_1\) の4次式と \(x_3\) の和が定数値と等しいことを表しています。
\(x_2, x_3\) の下限が0に設定されているのに対し、\(x_0, x_1\) は上下限が設定されていないことから、\(x_0, x_1\) は負の値も取り得ることがわかります。
Exampleの実行コマンド:
python -m naginterfaces.library.examples.opt.handle_solve_ssqp_exソースコード表示コマンド:
python -c "import inspect; from naginterfaces.library.examples.opt import handle_solve_ssqp_ex; print(''.join(inspect.getsourcelines(handle_solve_ssqp_ex)[0]))"出力結果例:
naginterfaces.library.opt.handle_solve_ssqp Python Example Results.
At the solution the objective function is 1.900124e+00.
マニュアル:
ソース:
#!/usr/bin/env python3
"``naginterfaces.library.opt.handle_solve_ssqp`` Python Example."
# nAG Copyright 2022.
# pylint: disable=invalid-name,too-many-arguments,too-many-locals,too-many-statements
import numpy as np
from naginterfaces.library import opt
def main():
"""
Example for :func:`naginterfaces.library.opt.handle_solve_ssqp`.
SQP solver for sparse NLP.
NLP example: Quadratic objective, linear constraint and two nonlinear
constraints.
For illustrative purposes, the quadratic objective is coded as a nonlinear
function to show the usage of objfun, objgrd user call-backs.
>>> main()
naginterfaces.library.opt.handle_solve_ssqp Python Example Results.
At the solution the objective function is 1.900124e+00.
"""
print(
'naginterfaces.library.opt.handle_solve_ssqp Python Example Results.'
)
# The nonlinear objective:
cb_objfun = lambda x, inform: (
(
(x[0]+x[1]+x[2])*(x[0]+x[1]+x[2]) + 3.0*x[2] + 5.0*x[3] +
np.cos(0.01*x[0]) - 1.0
),
0,
)
def cb_objgrd(x, fdx, inform):
"""The gradient of the objective."""
summ = 2.0*(x[0]+x[1]+x[2])
fdx[0] = summ - 0.01*np.sin(0.01*x[0])
fdx[1] = summ
fdx[2] = summ + 3.0
fdx[3] = 5.0
return 0
# The nonlinear constraints function:
cb_confun = lambda x, _ncnln, inform: (
[
x[0]*x[0] + x[1]*x[1] + x[2], x[1]*x[1]*x[1]*x[1] + x[3],
],
0,
)
def cb_congrd(x, gdx, inform):
"""The Jacobian of the nonlinear constraints."""
gdx[:] = [
2.0*x[0],
2.0*x[1],
1.0,
4.0*x[1]*x[1]*x[1],
1.0,
]
return inform
# The initial guess:
x = [1., 2., 3., 4.]
nvar = len(x)
nnzu = 2*nvar
# Create a handle for the problem:
handle = opt.handle_init(nvar)
# Define the bounds:
opt.handle_set_simplebounds(
handle,
bl=[-1.e20, -1.e20, 0., 0.], bu=[1.e20]*nvar,
)
# Define the non-linear objective:
opt.handle_set_nlnobj(
handle,
idxfd=[1, 2, 3, 4],
)
# Define the linear constraints:
opt.handle_set_linconstr(
handle,
bl=[0.], bu=[1.e20],
irowb=[1, 1],
icolb=[1, 2],
b=[2.0, 4.0],
)
nnzu += 2*1
# Define the nonlinear constraints:
opt.handle_set_nlnconstr(
handle,
bl=[2., 4.], bu=[2., 4.],
irowgd=[1, 1, 1, 2, 2], icolgd=[1, 2, 3, 2, 4],
)
nnzu += 2*2
opt.handle_opt_set(
handle, 'Verify Derivatives = Yes',
)
opt.handle_opt_set(
handle, 'Print Solution = No',
)
opt.handle_opt_set(
handle, 'Print Level = 0',
)
# Optionally, define variable x(4) to be linear
# Hinting which variables are linear in problems with many
# variables can speed up performance
opt.handle_set_property(handle, "Linear", 4)
# Solve the problem:
sn_soln = opt.handle_solve_ssqp(
handle, x,
objfun=cb_objfun, objgrd=cb_objgrd, confun=cb_confun, congrd=cb_congrd
)
print(
'At the solution the objective function is {:.6e}.'.format(
sn_soln.rinfo[0]
)
)
# Destroy the handle:
opt.handle_free(handle)
if __name__ == '__main__':
import doctest
import sys
sys.exit(
doctest.testmod(
None, verbose=True, report=False,
optionflags=doctest.ELLIPSIS | doctest.REPORT_NDIFF,
).failed
)