openMDAO: Does the use of ExecComp maximum() interfere with constraints not being affected by design variables?

Viewed 23

When running the optimization driver on a large model I recieve:

DerivativesWarning:Constraints or objectives [('max_current.current_constraint.current_constraint', inds=[0]), ('max_current.continuous_current_constraint.continuous_current_constraint', inds=[0])] cannot be impacted by the design variables of the problem.

I read the answer to a similar question posed here.

The values do change as the design variables change, and the two constraints are satisfied during the course of optimization.

I had assumed this was due to those components' ExecComp using a maximum(), as this is the only place in my model I use a maximum function, however when setting up a simple problem with a maximum() function in a similar manner I do not receive an error.

My model uses explicit components that are looped, there are connections in the bottom left of the N2 diagram and NLBGS is converging the whole model. I currently am thinking it is due to the use of only explicit components and the NLBGS instead of implicit components.

Thank you for any insight you can give in resolving this warning.

Below is a simple script using maximum() that does not report errors. (I was so sure that was it) As I create a minimum working example that gives the error in a similar way to my larger model I will upload it.

import openmdao.api as om

prob=om.Problem()
prob.driver = om.ScipyOptimizeDriver()
prob.driver.options['optimizer'] = 'SLSQP'
prob.driver.options['tol'] = 1e-6
prob.driver.options['maxiter'] = 80
prob.driver.options['disp'] = True

indeps = prob.model.add_subsystem('indeps', om.IndepVarComp())
indeps.add_output('x', val=2.0, units=None)
prob.model.promotes('indeps', outputs=['*'])

prob.model.add_subsystem('y_func_1',
                         om.ExecComp('y_func_1 = x'),
                         promotes_inputs=['x'],
                         promotes_outputs=['y_func_1'])
prob.model.add_subsystem('y_func_2',
                         om.ExecComp('y_func_2 = x**2'),
                         promotes_inputs=['x'],
                         promotes_outputs=['y_func_2'])
prob.model.add_subsystem('y_max',
                         om.ExecComp('y_max = maximum( y_func_1 , y_func_2 )'),
                         promotes_inputs=['y_func_1',
                                          'y_func_2'],
                         promotes_outputs=['y_max'])
prob.model.add_subsystem('y_check',
                         om.ExecComp('y_check = y_max - 1.1'),
                         promotes_inputs=['*'],
                         promotes_outputs=['*'])

prob.model.add_constraint('y_check', lower=0.0)

prob.model.add_design_var('x', lower=0.0, upper=2.0)
prob.model.add_objective('x')

prob.setup()
prob.run_driver()

print(prob.get_val('x'))
1 Answers

There is a problem with the maximum function in this context. Technically a maximum function is not differentiable; at least not when the index of which value is max is subject to change. If the maximum value is not subject to change, then it is differentiable... but you didn't need the max function anyway.

One correct, differentiable way to handle a max when doing gradient based things is to use a KS function. OpenMDAO provides the KSComp which implements it. There are other kinds of functions (like p-norm that you could use as well).

However, even though maximum is not technically differentiable ... you can sort-of/kind-of get away with it. At least, numpy (which ExecComp uses) lets you apply complex-step differentiation to the maximum function and it seems to give a non-zero derivative. So while its not technically correct, you can maybe get rid of it. At least, its not likely to be the core of your problem.

You mention using NLBGS, and that you have components which are looped. Your test case is purely feed forward though (here is the N2 from your test case). OpenMDAO N2 diagram of feed forward model. That is an important difference.

The problem here is with your derivatives, not with the maximum function. Since you have a nonlinear solver, you need to do something to get the derivatives right. In the example Sellar optimization, the model uses this line: prob.model.approx_totals(), which tells OpenMDAO to finite-difference across the whole model (including the nonlinear solver). This is simple and keeps the example compact. It also works regardless of whether your components define derivatives or not. It is however, slow and suffers from numerical difficulties. So use on "real" problems at your own risk.

If you don't include that (and your above example does not, so I assume your real problem does not either) then you're basically telling OpenMDAO that you want to use analytic derivatives (yay! they are so much more awesome). That means that you need to have a Linear solver to match your nonlinear one. For most problems that you start out with, you can simply put a DirectSolver right at the top of the model and it will all work out. For more advanced models, you need a more complex linear solver structure... but thats a different question.

Give this a try: prob.model.linear_solver = om.DirectSolver()

That should give you non-zero total derivatives regardless of whether you have coupling (loops) or not.

Related