1

I am implementing Bryson-Denham problem. The objective is:

$J=\frac{1}{2}\int_{0}^{1}u^2\left(t\right)dt$

and in the doc of Dymos, all explanation and examples state objective value as a scalar at loc=initial or loc=final. I could not find any example that use integral of some quantity over time as an objective function. Is this possible? How can I implement this?

FYI, Bryson-Denham problem is well-explained in this page: https://www.gpops2.com/Examples/Bryson-Denham.html

Konrad Rudolph
  • 530,221
  • 131
  • 937
  • 1,214
ehehd
  • 13
  • 3
  • Dymos considers any quantity integrated wrt time to be a state. Since the objective J is integrated on the same interval as the other states, we can just add it as a state itself. – Rob Falck Dec 16 '21 at 21:20
  • Rob, after I saw your comment and Justin's detailed explanation, I think I got a solid idea how Dymos works. Thank you a lot for your answer. – ehehd Dec 16 '21 at 22:17
  • Because this is such a classic problem I've added it to our examples in the documentation here: https://openmdao.github.io/dymos/examples/bryson_denham/bryson_denham.html – Rob Falck Dec 21 '21 at 21:20

1 Answers1

1

Dymos will integrate any state you give it. In this case, you need to add a state for J and then also compute a state rate for it --- J_dot.

import openmdao.api as om 
import dymos as dm

class BrysonDedhamODE(om.ExplicitComponent):


    def initialize(self):
        self.options.declare('num_nodes', types=int)

    def setup(self):
        nn = self.options['num_nodes']

        # static parameters
        self.add_input('x', shape=nn)
        self.add_input('v', shape=nn)
        self.add_input('u', shape=nn)
        self.add_input('J', shape=nn)



        # state rates
        self.add_output('x_dot', shape=nn, tags=['dymos.state_rate_source:x'])
        self.add_output('v_dot', shape=nn, tags=['dymos.state_rate_source:v'])
        self.add_output('J_dot', shape=nn, tags=['dymos.state_rate_source:J'])

        # Ask OpenMDAO to compute the partial derivatives using complex-step
        # with a partial coloring algorithm for improved performance, and use
        # a graph coloring algorithm to automatically detect the sparsity pattern.
        self.declare_coloring(wrt='*', method='cs')



    def compute(self, inputs, outputs):

        v, u, j = inputs["v"], inputs["u"], inputs["J"]

        outputs['x_dot'] = v        
        outputs['v_dot'] = u        
        outputs['J_dot'] = 0.5*u**2


p = om.Problem()

p.driver = om.pyOptSparseDriver()
p.driver.options['optimizer'] = 'SLSQP'
p.driver.declare_coloring()

traj = p.model.add_subsystem('traj', dm.Trajectory())

transcription = dm.Radau(num_segments=10, order=3, compressed=True)
phase0 = dm.Phase(ode_class=BrysonDedhamODE, transcription=transcription)

traj.add_phase('phase0', phase0)


phase0.set_time_options(fix_initial=True, fix_duration=True)

phase0.set_state_options("x", fix_initial=True, fix_final=True, lower=0, upper=2)
phase0.set_state_options("v", fix_initial=True, fix_final=True, lower=-2, upper=2)
phase0.set_state_options("J", fix_initial=False, fix_final=False,lower=-10, upper=10)

phase0.add_control('u', lower=-10, upper=10)
phase0.add_path_constraint('x', upper=1/9.)

phase0.add_objective('J', loc="final")

p.setup()

#initial conditions
p['traj.phase0.states:x'] = phase0.interp('x', [0,0])
p['traj.phase0.states:x'] = phase0.interp('x', [0,0])
p['traj.phase0.states:v'] = phase0.interp('v', [1,-1])
p['traj.phase0.t_duration'] = 1
p['traj.phase0.t_initial'] = 0

dm.run_problem(p, make_plots=True)

state history for v state history for x control history for u

Justin Gray
  • 5,605
  • 1
  • 11
  • 16
  • Justin's answer here is good, though you might want to specify `fix_initial=True` for J and set the initial value to zero. Also, you can specify `'u'` as the rate_source for `'v'`, and `'v'` as the rate_source for `'x'`. The ODE component only needs to compute `J_dot` in this case. – Rob Falck Dec 16 '21 at 21:22
  • Justin, thank you so much! This answer really helps a lot! – ehehd Dec 16 '21 at 22:15
  • Also thank you for additional comments, Rob. – ehehd Dec 16 '21 at 22:16