密な非線形計画問題の解法とアルゴリズムパラメータの設定

nAG数値計算ライブラリ > 最適化アルゴリズムExample集 > 密な非線形計画問題の解法とアルゴリズムパラメータの設定

非線形計画法-NLP（密）

このExampleは、Hock and Schittkowski問題71をNLP（非線形計画法）で解くことを目的としています。また、オプションのアルゴリズム・パラメータの取り扱い方を示しています。

目的関数:

タスク	式
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\)

追加の詳細説明:

コードでは、無限大の境界値を\(10^{20}\)と設定しています。これは”Infinite Bound Size”オプションで指定されています。
目的関数と制約条件の勾配(ヤコビアン)は、コールバック関数cb_objfunとcb_confun内で計算されています。
初期推定値は\(x = [1, 5, 5, 1]\)です。
最終的な目的関数値は約17.014017となります。

Exampleの実行コマンド:

python -m naginterfaces.library.examples.opt.nlp1_solve_ex

ソースコード表示コマンド:

python -c "import inspect; from naginterfaces.library.examples.opt import nlp1_solve_ex; print(''.join(inspect.getsourcelines(nlp1_solve_ex)[0]))"

出力結果例:

naginterfaces.library.opt.nlp1_solve Python Example Results.
Solve Hock and Schittkowski Problem 71.
Final objective value is 1.7014017e+01

マニュアル:

nlp1_solveのマニュアル

ソース:

#!/usr/bin/env python3
"``naginterfaces.library.opt.nlp1_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.nlp1_solve`.

    Dense NLP.

    Demonstrates handling optional algorithmic parameters.

    >>> main()
    naginterfaces.library.opt.nlp1_solve Python Example Results.
    Solve Hock and Schittkowski Problem 71.
    Final objective value is 1.7014017e+01
    """

    print(
        'naginterfaces.library.opt.nlp1_solve Python Example Results.'
    )
    print('Solve Hock and Schittkowski Problem 71.')

    def cb_confun(mode, needc, x, cjac, _nstate):
        """The nonlinear constraints."""
        c = np.zeros(len(needc))
        if needc[0] > 0:
            if mode in [0, 2]:
                c[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]:
                c[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 c, cjac

    def cb_objfun(mode, x, objgrd, _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:
            objgrd[:] = [
                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 objf, objgrd

    # Initialize the solver:
    comm = opt.nlp1_init('nlp1_solve')

    # The initial guess:
    x = [1., 5., 5., 1.]
    # The linear constraints:
    a = np.array([[1.]*len(x)])
    # The bounds:
    bl = [1., 1., 1., 1., -1.0E+25, -1.0E+25, 25.]
    bu = [5., 5., 5., 5., 20., 40., 1.0E+25]

    # To set algorithmic options:
    opt.nlp1_option_string('Infinite Bound Size = 1.0e20', comm)

    objf = opt.nlp1_solve(
        a, bl, bu, cb_objfun, 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
    )