Simple Optimization using Simultaneous Derivatives¶

Consider a set of points in 2-d space that are to be arranged along a circle such that the radius of the circle is maximized, subject to constraints. We start out with the points randomly distributed within a unit circle centered about the origin. The locations of our points are determined by the values of the x and y arrays defined in our problem.

from openmdao.api import Problem, IndepVarComp, ExecComp, ScipyOptimizeDriver
import numpy as np

# note: size must be an even number
SIZE = 10
p = Problem()

indeps = p.model.add_subsystem('indeps', IndepVarComp(), promotes_outputs=['*'])

# the following were randomly generated using np.random.random(10)*2-1 to randomly
# disperse them within a unit circle centered at the origin.
indeps.add_output('x', np.array([ 0.55994437, -0.95923447,  0.21798656, -0.02158783,  0.62183717,
0.04007379,  0.46044942, -0.10129622,  0.27720413, -0.37107886]))
indeps.add_output('y', np.array([ 0.52577864,  0.30894559,  0.8420792 ,  0.35039912, -0.67290778,
-0.86236787, -0.97500023,  0.47739414,  0.51174103,  0.10052582]))

p.model.add_subsystem('r_con', ExecComp('g=x**2 + y**2 - r',
g=np.ones(SIZE), x=np.ones(SIZE), y=np.ones(SIZE)))

thetas = np.linspace(0, np.pi/4, SIZE)
p.model.add_subsystem('theta_con', ExecComp('g=arctan(y/x) - theta',
g=np.ones(SIZE), x=np.ones(SIZE), y=np.ones(SIZE),
theta=thetas))
p.model.add_subsystem('delta_theta_con', ExecComp('g = arctan(y/x)[::2]-arctan(y/x)[1::2]',
g=np.ones(SIZE//2), x=np.ones(SIZE),
y=np.ones(SIZE)))

thetas = np.linspace(0, np.pi/4, SIZE)

p.model.add_subsystem('l_conx', ExecComp('g=x-1', g=np.ones(SIZE), x=np.ones(SIZE)))

p.model.connect('r', ('circle.r', 'r_con.r'))
p.model.connect('x', ['r_con.x', 'theta_con.x', 'delta_theta_con.x'])

p.model.connect('x', 'l_conx.x')

p.model.connect('y', ['r_con.y', 'theta_con.y', 'delta_theta_con.y'])

p.driver = ScipyOptimizeDriver()
p.driver.options['optimizer'] = 'SLSQP'
p.driver.options['disp'] = False

# nonlinear constraints

IND = np.arange(SIZE, dtype=int)
ODD_IND = IND[0::2]  # all odd indices
p.model.add_constraint('theta_con.g', lower=-1e-5, upper=1e-5, indices=ODD_IND)

# this constrains x[0] to be 1 (see definition of l_conx)
p.model.add_constraint('l_conx.g', equals=0, linear=False, indices=[0,])

# linear constraint
p.model.add_constraint('y', equals=0, indices=[0,], linear=True)

# setup coloring
p.driver.set_simul_deriv_color(
(
{
'y': [0, 1, 0, 1, 0, 1, 0, 1, 0, 1],
'x': [0, 1, 0, 1, 0, 1, 0, 1, 0, 1]
},
{
'delta_theta_con.g': {
'y': {
0: ([0, 1, 2, 3, 4], [0, 2, 4, 6, 8]),
1: ([0, 1, 2, 3, 4], [1, 3, 5, 7, 9])},
'x': {
0: ([0, 1, 2, 3, 4], [0, 2, 4, 6, 8]),
1: ([0, 1, 2, 3, 4], [1, 3, 5, 7, 9])}},
'r_con.g': {
'y': {
0: ([0, 2, 4, 6, 8], [0, 2, 4, 6, 8]),
1: ([1, 3, 5, 7, 9], [1, 3, 5, 7, 9])},
'x': {
0: ([0, 2, 4, 6, 8], [0, 2, 4, 6, 8]),
1: ([1, 3, 5, 7, 9], [1, 3, 5, 7, 9])}},
'l_conx.g': {
'x': {
0: ([0], [0])}},
'theta_con.g': {
'y': {
0: ([0, 1, 2, 3, 4], [0, 2, 4, 6, 8])},
'x': {
0: ([0, 1, 2, 3, 4], [0, 2, 4, 6, 8])
}
}
}
)
)

p.setup(mode='fwd')
p.run_driver()

print(p['circle.area'])

[ 3.14159265]

Total derivatives with respect to x and y will be solved for simultaneously based on the color of the points shown below. Note that we have two colors and our x and y design variables are of size 10. We have a third design variable, r, that is size 1. This means that if we don’t solve for derivatives simultaneously, we must perform 21 linear solves (10 + 10 + 1) to solve for total derivatives with respect to all of our design variables. But, since each of our x and y design variables have only two colors, we can solve for all of our total derivatives using only 5 linear solves (2 + 2 + 1). This means that using simultaneous derivatives saves us 16 linear solves each time we compute our total derivatives.

Here’s our problem at the start of the optimization:

After the optimization, all of our points lie along the unit circle. The final radius is 1.0 (which corresponds to an area of PI) because we constrained our x[0] value to be 1.0.