疎な非線形最適化問題の逐次二次計画法による解法

nAG数値計算ライブラリ > 最適化アルゴリズムExample集 > 疎な非線形最適化問題の逐次二次計画法による解法

逐次二次計画法-SQP（疎）

このExampleは、スパースな非線形計画問題(NLP)を逐次二次計画法(SQP)で解いています。2次の目的関数、線形制約、2つの非線形制約を含む問題を扱っています。目的関数は、ユーザコールバックの使用方法を示すために非線形関数としてコード化されています。

目的関数:

タスク	式
minimize	$f (x) = (x_{0} + x_{1} + x_{2})^{2} + 3 x_{2} + 5 x_{3} + \cos (0.01 x_{0}) - 1$

決定変数:

変数	範囲
$x_{0}$	$(- \infty, \infty)$
$x_{1}$	$(- \infty, \infty)$
$x_{2}$	$[0, \infty)$
$x_{3}$	$[0, \infty)$

制約条件:

制約	式
線形制約	$2 x_{0} + 4 x_{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.01 x_{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.

マニュアル:

handle_solve_ssqpのマニュアル

ソース:

#!/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
    )