1.3. Example 4
An example written as a Python notebook (.py) with minimal explanation in Markdown
We start this script by defining the optimization problem.
'''Example 4'''
import numpy as np
from modopt.api import Problem
class X4(Problem):
def initialize(self, ):
# Name your problem
self.problem_name = 'x^4'
def setup(self):
# Add design variables of your problem
self.add_design_variables('x',
shape=(2, ),
vals=np.array([.3, .3]))
self.add_objective('f',)
def setup_derivatives(self):
# Declare objective gradient and its shape
self.declare_objective_gradient(wrt='x', )
self.declare_objective_hessian(of='x', wrt='x')
# Compute the value of the objective, gradient and Hessian
# with the given design variable values
def compute_objective(self, dvs, obj):
obj['f'] = np.sum(dvs['x']**4)
def compute_objective_gradient(self, dvs, grad):
grad['x'] = 4 * dvs['x']**3
def compute_objective_hessian(self, dvs, hess):
hess['x', 'x'] = 12 * np.diag(dvs['x']**2)
We will now build the steepest descent optimization algorithm using the optimizer class.
import numpy as np
import time
from modopt.api import Optimizer
class SteepestDescent(Optimizer):
def initialize(self):
# Name your algorithm
self.solver_name = 'steepest_descent'
self.obj = self.problem._compute_objective
self.grad = self.problem._compute_objective_gradient
self.options.declare('max_itr', default=1000, types=int)
self.options.declare('opt_tol', types=float)
# Specify format of outputs available from your optimizer after each iteration
self.default_outputs_format = {
'itr': int,
'obj': float,
# for arrays from each iteration, shapes need to be declared
'x': (float, (self.problem.nx, )),
'opt': float,
'time': float,
}
# Enable user to specify, as a list, which among the available outputs
# need to be stored in memory and written to output files
self.options.declare('outputs',
types=list,
default=['itr', 'obj', 'x', 'opt', 'time'])
def solve(self):
nx = self.problem.nx
x = self.problem.x.get_data()
opt_tol = self.options['opt_tol']
max_itr = self.options['max_itr']
obj = self.obj
grad = self.grad
start_time = time.time()
# Setting intial values for initial iterates
x_k = x * 1.
f_k = obj(x_k)
g_k = grad(x_k)
# Iteration counter
itr = 0
# Optimality
opt = np.linalg.norm(g_k)
# Initializing outputs
self.update_outputs(itr=0,
x=x_k,
obj=f_k,
opt=opt,
time=time.time() - start_time)
while (opt > opt_tol and itr < max_itr):
itr_start = time.time()
itr += 1
# ALGORITHM STARTS HERE
# >>>>>>>>>>>>>>>>>>>>>
p_k = -g_k
x_k += p_k
f_k = obj(x_k)
g_k = grad(x_k)
opt = np.linalg.norm(g_k)
# <<<<<<<<<<<<<<<<<<<
# ALGORITHM ENDS HERE
# Append arrays inside outputs dict with new values from the current iteration
self.update_outputs(itr=itr,
x=x_k,
obj=f_k,
opt=opt,
time=time.time() - start_time)
# Run post-processing for the Optimizer() base class
self.run_post_processing()
end_time = time.time()
self.total_time = end_time - start_time
Now set up your optimizer with the problem defined above to solve it.
# Set your optimality tolerance
opt_tol = 1E-8
# Set maximum optimizer iteration limit
max_itr = 100
prob = X4()
from modopt.optimization_algorithms import Newton, QuasiNewton, SQP
# Set up your optimizer with your problem and pass in optimizer parameters
optimizer = SteepestDescent(prob,
opt_tol=opt_tol,
max_itr=max_itr,
outputs=['itr', 'obj', 'x', 'opt', 'time'])
optimizer = Newton(prob, opt_tol=opt_tol)
optimizer = QuasiNewton(prob, opt_tol=opt_tol)
# Check first derivatives at the initial guess, if needed
optimizer.check_first_derivatives(prob.x.get_data())
# Solve your optimization problem
optimizer.solve()
# Print results of optimization (summary_table contains information from each iteration)
optimizer.print_results(summary_table=True, compact_print=True)
Setting objective name as "f".
Directory x^4_outputs already exists
Directory x^4_outputs already exists
----------------------------------------------------------------------------
Derivative type | Calc norm | FD norm | Abs error norm | Rel error norm
----------------------------------------------------------------------------
Gradient | 1.5274e-01 | 1.5274e-01 | 7.6367e-07 | 7.0710e-06
----------------------------------------------------------------------------
===============================
ModOpt final iteration summary:
===============================
Problem : x^4
Solver : bfgs
itr : 20
obj : 1.6785177544170953e-12
opt : 4.9601744115117986e-09
time : 0.403872013092041
num_f_evals : 21
num_g_evals : 21
step : 1.0
=========================================
==================================================================
modOpt summary table:
==================================================================
itr obj opt time num_f_evals num_g_evals step
0 1.62E-02 1.53E-01 5.67E-05 1 1 0.00E+00
2 1.11E-03 2.05E-02 3.28E-02 3 3 1.00E+00
4 1.11E-04 3.64E-03 6.53E-02 5 5 1.00E+00
6 1.16E-05 6.69E-04 8.97E-02 7 7 1.00E+00
8 1.22E-06 1.24E-04 1.18E-01 9 9 1.00E+00
11 4.18E-08 9.84E-06 1.90E-01 12 12 1.00E+00
13 4.41E-09 1.82E-06 2.28E-01 14 14 1.00E+00
15 4.65E-10 3.37E-07 2.96E-01 16 16 1.00E+00
17 4.90E-11 6.23E-08 3.43E-01 18 18 1.00E+00
20 1.68E-12 4.96E-09 4.04E-01 21 21 1.00E+00
==================================================================
Print to console to see any of the outputs that were declared. Since the arrays are long, here we only print the last entry and verify it with the print_results() above.
print('\n')
print('Optimizer data')
print('num_iterations:', optimizer.outputs['itr'][-1])
print('optimized_dvs:', optimizer.outputs['x'][-1])
print('optimization_time:', optimizer.outputs['time'][-1])
print('optimized_obj:', optimizer.outputs['obj'][-1])
print('final_optimality:', optimizer.outputs['opt'][-1])
print('\n')
print('Final problem data')
print('optimized_dvs:', prob.x.get_data())
print('optimized_obj:', prob.obj['f'])
Optimizer data
num_iterations: 20
optimized_dvs: [0.00095711 0.00095716]
optimization_time: 0.403872013092041
optimized_obj: 1.6785177544170953e-12
final_optimality: 4.9601744115117986e-09
Final problem data
optimized_dvs: [0.00095711 0.00095716]
optimized_obj: 1.6785177544170953e-12