Source code for naginterfaces.library.examples.opt.handle_solve_ipopt_ex

#!/usr/bin/env python3
"``naginterfaces.library.opt.handle_solve_ipopt`` Python Example."

# NAG Copyright 2018-2020.

# pylint: disable=invalid-name,too-many-arguments,too-many-locals,too-many-statements

from math import sqrt

import numpy as np

from naginterfaces.library import opt

[docs]def main(): """ Example for :func:`naginterfaces.library.opt.handle_solve_ipopt`. Interior-point solver for sparse NLP. >>> main() naginterfaces.library.opt.handle_solve_ipopt Python Example Results. Solving a problem based on Hock and Schittkowski Problem 73. Solving with a nonlinear objective. At the solution the objective function is 2.9894378e+01. """ print( 'naginterfaces.library.opt.handle_solve_ipopt Python Example Results.' ) print('Solving a problem based on Hock and Schittkowski Problem 73.') print('Solving with a nonlinear objective.') # The 'non'linear objective: cb_objfun = lambda x, inform: ( 24.55*x[0] + 26.75*x[1] + 39.*x[2] + 40.5*x[3], inform, ) def cb_objgrd(x, fdx, inform): fdx[:] = [24.55, 26.75, 39., 40.5] return inform cb_confun = lambda x, _ncnln, inform: ( [ 12.0*x[0] + 11.9*x[1] + 41.8*x[2] + 52.1*x[3] - 1.645 * sqrt( 0.28*x[0]**2 + 0.19*x[1]**2 + 20.5*x[2]**2 + 0.62*x[3]**2 ) ], inform, ) def cb_congrd(x, gdx, inform): """The Jacobian of the nonlinear constraints.""" tmp = sqrt( 0.62*x[3]**2 + 20.5*x[2]**2 + 0.19*x[1]**2 + 0.28*x[0]**2 ) gdx[:] = [ (12.0*tmp-0.4606*x[0])/tmp, (11.9*tmp-0.31255*x[1])/tmp, (41.8*tmp-33.7225*x[2])/tmp, (52.1*tmp-1.0199*x[3])/tmp, ] return inform def cb_hess(x, idf, sigma, lamda, hx, inform): """Hessian function.""" nnzh = hx.size if idf == 0: # Objective is linear. for i in range(nnzh): hx[i] = 0.0 return inform tmp = sqrt( 0.62*x[3]**2 + 20.5*x[2]**2 + 0.19*x[1]**2 + 0.28*x[0]**2 ) tmp = tmp*( x[3]**2 + 33.064516129032258064*x[2]**2 + 0.30645161290322580645*x[1]**2 + 0.45161290322580645161*x[0]**2 ) hx[0] = ( -0.4606*x[3]**2 - 15.229516129032258064*x[2]**2 - 0.14115161290322580645*x[1]*x[1])/tmp hx[1] = (0.14115161290322580645*x[0]*x[1])/tmp hx[2] = (15.229516129032258064*x[0]*x[2])/tmp hx[3] = (0.4606*x[0]*x[3])/tmp hx[4] = ( -0.31255*x[3]**2 - 10.334314516129032258*x[2]**2 - 0.14115161290322580645*x[0]**2 )/tmp hx[5] = (10.334314516129032258*x[1]*x[2])/tmp hx[6] = (0.31255*x[1]*x[3])/tmp hx[7] = ( -33.7225*x[3]**2 - 10.334314516129032258*x[1]**2 - 15.229516129032258065*x[0]**2 )/tmp hx[8] = (33.7225*x[2]*x[3])/tmp hx[9] = ( -33.7225*x[2]**2 - 0.31255*x[1]**2 - 0.4606*x[0]**2 )/tmp if idf == -1: for i in range(nnzh): hx[i] = hx[i] * lamda[0] else: assert idf == 1 return inform # The initial guess: x = [1., 1., 1., 1.] 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=[0.]*nvar, 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=[5., 1.], bu=[1.e20, 1.], irowb=[1, 1, 1, 1, 2, 2, 2, 2], icolb=[1, 2, 3, 4, 1, 2, 3, 4], b=[2.3, 5.6, 11.1, 1.3, 1.0, 1.0, 1.0, 1.0], ) nnzu += 2*2 # Define the nonlinear constraints: opt.handle_set_nlnconstr( handle, bl=[21.], bu=[1.e20], irowgd=[1, 1, 1, 1], icolgd=[1, 2, 3, 4], ) nnzu += 2*1 # Define the structure of the Hessian, dense upper triangle: irowh = np.empty(10, dtype=int) icolh = np.empty(10, dtype=int) idx = 0 for i in range(1, nvar + 1): for j in range(i, nvar + 1): icolh[idx] = j irowh[idx] = i idx += 1 # For the objective function: opt.handle_set_nlnhess( handle, idf=0, irowh=irowh, icolh=icolh, ) # For the constraint function: opt.handle_set_nlnhess( handle, idf=1, irowh=irowh, icolh=icolh, ) opt.handle_opt_set( handle, 'Stop Tolerance 1 = 2.5e-8', ) opt.handle_opt_set( handle, 'Print Level = 0', ) # Solve the problem: ip_soln = opt.handle_solve_ipopt( handle, x, objfun=cb_objfun, objgrd=cb_objgrd, confun=cb_confun, congrd=cb_congrd, hess=cb_hess, ) print( 'At the solution the objective function is {:.7e}.'.format( ip_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 )