Cannonball Implicit Duration#
+This example demonstrates determining the duration of a phase using a +nonlinear solver to satisfy an endpoint condition. As in the multi-phase +cannonball problem, we will simultaneously find the optimal design for +the cannonball (its radius) and the optimal flight profile (its launch +angle). However, here we will use a single phase that terminates at the +condition
+The cannonball sizing component#
+This component computes the area (needed to compute drag) and mass (needed to compute energy) of a cannonball of a given radius and density.
+This component sits upstream of the trajectory model and feeds its outputs to the trajectory as parameters.
+import numpy as np
+from scipy.interpolate import interp1d
+import matplotlib.pyplot as plt
+
+import openmdao.api as om
+
+import dymos as dm
+from dymos.models.atmosphere.atmos_1976 import USatm1976Data
+
+#############################################
+# Component for the design part of the model
+#############################################
+class CannonballSizeComp(om.ExplicitComponent):
+ """
+ Compute the area and mass of a cannonball with a given radius and density.
+
+ Notes
+ -----
+ This component is not vectorized with 'num_nodes' as is the usual way
+ with Dymos, but is instead intended to compute a scalar mass and reference
+ area from scalar radius and density inputs. This component does not reside
+ in the ODE but instead its outputs are connected to the trajectory via
+ input design parameters.
+ """
+ def setup(self):
+ self.add_input(name='radius', val=1.0, desc='cannonball radius', units='m')
+ self.add_input(name='dens', val=7870., desc='cannonball density', units='kg/m**3')
+
+ self.add_output(name='mass', shape=(1,), desc='cannonball mass', units='kg')
+ self.add_output(name='S', shape=(1,), desc='aerodynamic reference area', units='m**2')
+
+ self.declare_partials(of='mass', wrt='dens')
+ self.declare_partials(of='mass', wrt='radius')
+
+ self.declare_partials(of='S', wrt='radius')
+
+ def compute(self, inputs, outputs):
+ radius = inputs['radius']
+ dens = inputs['dens']
+
+ outputs['mass'] = (4/3.) * dens * np.pi * radius ** 3
+ outputs['S'] = np.pi * radius ** 2
+
+ def compute_partials(self, inputs, partials):
+ radius = inputs['radius']
+ dens = inputs['dens']
+
+ partials['mass', 'dens'] = (4/3.) * np.pi * radius ** 3
+ partials['mass', 'radius'] = 4. * dens * np.pi * radius ** 2
+
+ partials['S', 'radius'] = 2 * np.pi * radius
+
The cannonball ODE component#
+This component computes the state rates and the kinetic energy of the cannonball.
+By calling the declare_coloring
method wrt all inputs and using method 'cs'
, we’re telling OpenMDAO to automatically determine the sparsity pattern of the outputs with respect to the inputs, and to automatically compute those outputs using complex-step approximation.
class CannonballODE(om.ExplicitComponent):
+ """
+ Cannonball ODE assuming flat earth and accounting for air resistance
+ """
+
+ def initialize(self):
+ self.options.declare('num_nodes', types=int)
+
+ def setup(self):
+ nn = self.options['num_nodes']
+
+ # static parameters
+ self.add_input('m', units='kg')
+ self.add_input('S', units='m**2')
+ # 0.5 good assumption for a sphere
+ self.add_input('CD', 0.5)
+
+ # time varying inputs
+ self.add_input('h', units='m', shape=nn)
+ self.add_input('v', units='m/s', shape=nn)
+ self.add_input('gam', units='rad', shape=nn)
+
+ # state rates
+ self.add_output('v_dot', shape=nn, units='m/s**2', tags=['dymos.state_rate_source:v'])
+ self.add_output('gam_dot', shape=nn, units='rad/s', tags=['dymos.state_rate_source:gam'])
+ self.add_output('h_dot', shape=nn, units='m/s', tags=['dymos.state_rate_source:h'])
+ self.add_output('r_dot', shape=nn, units='m/s', tags=['dymos.state_rate_source:r'])
+ self.add_output('ke', shape=nn, units='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')
+
+ alt_data = USatm1976Data.alt * om.unit_conversion('ft', 'm')[0]
+ rho_data = USatm1976Data.rho * om.unit_conversion('slug/ft**3', 'kg/m**3')[0]
+ self.rho_interp = interp1d(np.array(alt_data, dtype=complex),
+ np.array(rho_data, dtype=complex),
+ kind='linear')
+
+ def compute(self, inputs, outputs):
+
+ gam = inputs['gam']
+ v = inputs['v']
+ h = inputs['h']
+ m = inputs['m']
+ S = inputs['S']
+ CD = inputs['CD']
+
+ GRAVITY = 9.80665 # m/s**2
+
+ # handle complex-step gracefully from the interpolant
+ if np.iscomplexobj(h):
+ rho = self.rho_interp(inputs['h'])
+ else:
+ rho = self.rho_interp(inputs['h']).real
+
+ q = 0.5*rho*inputs['v']**2
+ qS = q * S
+ D = qS * CD
+ cgam = np.cos(gam)
+ sgam = np.sin(gam)
+ outputs['v_dot'] = - D/m-GRAVITY*sgam
+ outputs['gam_dot'] = -(GRAVITY/v)*cgam
+ outputs['h_dot'] = v*sgam
+ outputs['r_dot'] = v*cgam
+ outputs['ke'] = 0.5*m*v**2
+
Generating an initial guess#
+Since we will be using a nonlinear solver to converge the dynamics and the duration, a simple linear initial guess will not be sufficient. The initial guess we will use is generated from the kinematics equations that assume no atmospheric drag. We compute the state values for some duration and initial flight path angle
+def initial_guess(t_dur, gam0=np.pi/3):
+ t = np.linspace(0, t_dur, num=100)
+ g = -9.80665
+ v0 = -0.5*g*t_dur/np.sin(gam0)
+ r = v0*np.cos(gam0)*t
+ h = v0*np.sin(gam0)*t + 0.5*g*t**2
+ v = np.sqrt(v0*np.cos(gam0)**2 + (v0*np.sin(gam0) + g*t)**2)
+ gam = np.arctan2(v0*np.sin(gam0) + g*t, v0*np.cos(gam0)**2)
+
+ guess = {'t': t, 'r': r, 'h': h, 'v': v, 'gam': gam}
+
+ return guess
+
Building and running the problem#
+The following code defines the components for the physical +cannonball calculations and ODE problem, sets up trajectory using a +single phase, and specifies the stopping condition for the phase. +The initial flight path angle is free, since 45 degrees is not +necessarily optimal once air resistance is taken into account.
+p = om.Problem(model=om.Group())
+
+p.driver = om.pyOptSparseDriver()
+p.driver.options['optimizer'] = 'IPOPT'
+p.driver.declare_coloring()
+
+p.driver.opt_settings['derivative_test'] = 'first-order'
+p.driver.opt_settings['mu_strategy'] = 'monotone'
+p.driver.opt_settings['alpha_for_y'] = 'safer-min-dual-infeas'
+p.driver.opt_settings['bound_mult_init_method'] = 'mu-based'
+p.driver.opt_settings['mu_init'] = 0.01
+p.driver.opt_settings['nlp_scaling_method'] = 'gradient-based'
+
+p.set_solver_print(level=0, depth=99)
+
+p.model.add_subsystem('size_comp', CannonballSizeComp(),
+ promotes_inputs=['radius', 'dens'])
+p.model.set_input_defaults('dens', val=7.87, units='g/cm**3')
+p.model.add_design_var('radius', lower=0.01, upper=0.10,
+ ref0=0.01, ref=0.10, units='m')
+
+traj = p.model.add_subsystem('traj', dm.Trajectory())
+
+transcription = dm.Radau(num_segments=25, order=3, compressed=False)
+phase = dm.Phase(ode_class=CannonballODE, transcription=transcription)
+
+phase = traj.add_phase('phase', phase)
+
+# Duration is bounded to be greater than 1 to ensure the nonlinear solve
+# does not converge to the initial point.
+phase.set_time_options(fix_initial=True, duration_bounds=(1, 100), units='s')
+
+# All initial states except flight path angle are fixed
+# The output of the ODE which provides the rate source for each state
+# is obtained from the tags used on those outputs in the ODE.
+# The units of the states are automatically inferred by multiplying the units
+# of those rates by the time units.
+phase.add_state('r', fix_initial=True, solve_segments='forward')
+phase.add_state('h', fix_initial=True, solve_segments='forward')
+phase.add_state('gam', fix_initial=False, solve_segments='forward')
+phase.add_state('v', fix_initial=False, solve_segments='forward')
+
+phase.add_parameter('S', units='m**2', static_target=True)
+phase.add_parameter('m', units='kg', static_target=True)
+phase.add_parameter('CD', val=0.5, opt=False, static_target=True)
+
+phase.add_boundary_constraint('ke', loc='initial',
+ upper=400000, lower=0, ref=100000)
+
+# A duration balance is added setting altitude to zero.
+# A nonlinear solver is used to find the duration of required to satisfy the condition.
+phase.set_duration_balance('h', val=0.0)
+
+phase.add_objective('r', loc='final', scaler=-1.0)
+
+p.model.connect('size_comp.mass', 'traj.phase.parameters:m')
+p.model.connect('size_comp.S', 'traj.phase.parameters:S')
+
+# Finish Problem Setup
+p.setup()
+
+#############################################
+# Set constants and initial guesses
+#############################################
+p.set_val('radius', 0.05, units='m')
+p.set_val('dens', 7.87, units='g/cm**3')
+
+p.set_val('traj.phase.parameters:CD', 0.5)
+
+guess = initial_guess(t_dur=30.0)
+
+p.set_val('traj.phase.t_initial', 0.0)
+p.set_val('traj.phase.t_duration', 30.0)
+
+p.set_val('traj.phase.states:r', phase.interp('r', ys=guess['r'], xs=guess['t']))
+p.set_val('traj.phase.states:h', phase.interp('h', ys=guess['h'], xs=guess['t']))
+p.set_val('traj.phase.states:v', phase.interp('v', ys=guess['v'], xs=guess['t']))
+p.set_val('traj.phase.states:gam', phase.interp('gam', ys=guess['gam'], xs=guess['t']), units='rad')
+
+#####################################################
+# Run the optimization and final explicit simulation
+#####################################################
+dm.run_problem(p)
+
--- Constraint Report [traj] ---
+ --- phase ---
+ [initial] 0.0000e+00 <= ke <= 4.0000e+05 [J]
+
+Model viewer data has already been recorded for Driver.
+
Model viewer data has already been recorded for Driver.
+
Coloring for 'traj.phases.phase.rhs_all' (class CannonballODE)
+
+Jacobian shape: (500, 303) ( 0.92% nonzero)
+FWD solves: 6 REV solves: 0
+Total colors vs. total size: 6 vs 303 (98.0% improvement)
+
+Sparsity computed using tolerance: 1e-25
+Time to compute sparsity: 0.131124 sec.
+Time to compute coloring: 0.050527 sec.
+Memory to compute coloring: 0.000000 MB.
+
Full total jacobian was computed 3 times, taking 0.457985 seconds.
+Total jacobian shape: (98, 100)
+
+
+Jacobian shape: (98, 100) (21.30% nonzero)
+FWD solves: 24 REV solves: 0
+Total colors vs. total size: 24 vs 100 (76.0% improvement)
+
+Sparsity computed using tolerance: 1e-25
+Time to compute sparsity: 0.457985 sec.
+Time to compute coloring: 0.109271 sec.
+Memory to compute coloring: 1.652344 MB.
+
Optimization Problem -- Optimization using pyOpt_sparse
+================================================================================
+ Objective Function: _objfunc
+
+ Solution:
+--------------------------------------------------------------------------------
+ Total Time: 26.7678
+ User Objective Time : 23.0678
+ User Sensitivity Time : 3.2815
+ Interface Time : 0.2269
+ Opt Solver Time: 0.1916
+ Calls to Objective Function : 233
+ Calls to Sens Function : 54
+
+
+ Objectives
+ Index Name Value
+ 0 traj.phases.phase.indep_states.states:r -3.182962E+03
+
+ Variables (c - continuous, i - integer, d - discrete)
+ Index Name Type Lower Bound Value Upper Bound Status
+ 0 radius_0 c 0.000000E+00 3.540939E-01 1.000000E+00
+ 1 traj.phase.t_duration_0 c 1.000000E+00 1.000000E+02 1.000000E+02 u
+ 2 traj.phases.phase.indep_states.states:r_0 c -1.000000E+21 4.358669E+02 1.000000E+21
+ 3 traj.phases.phase.indep_states.states:r_1 c -1.000000E+21 7.590741E+02 1.000000E+21
+ 4 traj.phases.phase.indep_states.states:r_2 c -1.000000E+21 1.017857E+03 1.000000E+21
+ 5 traj.phases.phase.indep_states.states:r_3 c -1.000000E+21 1.234831E+03 1.000000E+21
+ 6 traj.phases.phase.indep_states.states:r_4 c -1.000000E+21 1.422453E+03 1.000000E+21
+ 7 traj.phases.phase.indep_states.states:r_5 c -1.000000E+21 1.588316E+03 1.000000E+21
+ 8 traj.phases.phase.indep_states.states:r_6 c -1.000000E+21 1.737381E+03 1.000000E+21
+ 9 traj.phases.phase.indep_states.states:r_7 c -1.000000E+21 1.873054E+03 1.000000E+21
+ 10 traj.phases.phase.indep_states.states:r_8 c -1.000000E+21 1.997754E+03 1.000000E+21
+ 11 traj.phases.phase.indep_states.states:r_9 c -1.000000E+21 2.113240E+03 1.000000E+21
+ 12 traj.phases.phase.indep_states.states:r_10 c -1.000000E+21 2.220813E+03 1.000000E+21
+ 13 traj.phases.phase.indep_states.states:r_11 c -1.000000E+21 2.321446E+03 1.000000E+21
+ 14 traj.phases.phase.indep_states.states:r_12 c -1.000000E+21 2.415862E+03 1.000000E+21
+ 15 traj.phases.phase.indep_states.states:r_13 c -1.000000E+21 2.504609E+03 1.000000E+21
+ 16 traj.phases.phase.indep_states.states:r_14 c -1.000000E+21 2.588097E+03 1.000000E+21
+ 17 traj.phases.phase.indep_states.states:r_15 c -1.000000E+21 2.666639E+03 1.000000E+21
+ 18 traj.phases.phase.indep_states.states:r_16 c -1.000000E+21 2.740483E+03 1.000000E+21
+ 19 traj.phases.phase.indep_states.states:r_17 c -1.000000E+21 2.809828E+03 1.000000E+21
+ 20 traj.phases.phase.indep_states.states:r_18 c -1.000000E+21 2.874850E+03 1.000000E+21
+ 21 traj.phases.phase.indep_states.states:r_19 c -1.000000E+21 2.935708E+03 1.000000E+21
+ 22 traj.phases.phase.indep_states.states:r_20 c -1.000000E+21 2.992555E+03 1.000000E+21
+ 23 traj.phases.phase.indep_states.states:r_21 c -1.000000E+21 3.045544E+03 1.000000E+21
+ 24 traj.phases.phase.indep_states.states:r_22 c -1.000000E+21 3.094833E+03 1.000000E+21
+ 25 traj.phases.phase.indep_states.states:r_23 c -1.000000E+21 3.140584E+03 1.000000E+21
+ 26 traj.phases.phase.indep_states.states:h_0 c -1.000000E+21 2.675246E+02 1.000000E+21
+ 27 traj.phases.phase.indep_states.states:h_1 c -1.000000E+21 4.559862E+02 1.000000E+21
+ 28 traj.phases.phase.indep_states.states:h_2 c -1.000000E+21 5.966683E+02 1.000000E+21
+ 29 traj.phases.phase.indep_states.states:h_3 c -1.000000E+21 7.042126E+02 1.000000E+21
+ 30 traj.phases.phase.indep_states.states:h_4 c -1.000000E+21 7.866604E+02 1.000000E+21
+ 31 traj.phases.phase.indep_states.states:h_5 c -1.000000E+21 8.488963E+02 1.000000E+21
+ 32 traj.phases.phase.indep_states.states:h_6 c -1.000000E+21 8.941007E+02 1.000000E+21
+ 33 traj.phases.phase.indep_states.states:h_7 c -1.000000E+21 9.244542E+02 1.000000E+21
+ 34 traj.phases.phase.indep_states.states:h_8 c -1.000000E+21 9.415152E+02 1.000000E+21
+ 35 traj.phases.phase.indep_states.states:h_9 c -1.000000E+21 9.464412E+02 1.000000E+21
+ 36 traj.phases.phase.indep_states.states:h_10 c -1.000000E+21 9.401270E+02 1.000000E+21
+ 37 traj.phases.phase.indep_states.states:h_11 c -1.000000E+21 9.232973E+02 1.000000E+21
+ 38 traj.phases.phase.indep_states.states:h_12 c -1.000000E+21 8.965701E+02 1.000000E+21
+ 39 traj.phases.phase.indep_states.states:h_13 c -1.000000E+21 8.605033E+02 1.000000E+21
+ 40 traj.phases.phase.indep_states.states:h_14 c -1.000000E+21 8.156262E+02 1.000000E+21
+ 41 traj.phases.phase.indep_states.states:h_15 c -1.000000E+21 7.624607E+02 1.000000E+21
+ 42 traj.phases.phase.indep_states.states:h_16 c -1.000000E+21 7.015341E+02 1.000000E+21
+ 43 traj.phases.phase.indep_states.states:h_17 c -1.000000E+21 6.333833E+02 1.000000E+21
+ 44 traj.phases.phase.indep_states.states:h_18 c -1.000000E+21 5.585557E+02 1.000000E+21
+ 45 traj.phases.phase.indep_states.states:h_19 c -1.000000E+21 4.776046E+02 1.000000E+21
+ 46 traj.phases.phase.indep_states.states:h_20 c -1.000000E+21 3.910839E+02 1.000000E+21
+ 47 traj.phases.phase.indep_states.states:h_21 c -1.000000E+21 2.995404E+02 1.000000E+21
+ 48 traj.phases.phase.indep_states.states:h_22 c -1.000000E+21 2.035072E+02 1.000000E+21
+ 49 traj.phases.phase.indep_states.states:h_23 c -1.000000E+21 1.034977E+02 1.000000E+21
+ 50 traj.phases.phase.indep_states.states:gam_0 c -1.000000E+21 5.588424E-01 1.000000E+21
+ 51 traj.phases.phase.indep_states.states:gam_1 c -1.000000E+21 5.398186E-01 1.000000E+21
+ 52 traj.phases.phase.indep_states.states:gam_2 c -1.000000E+21 5.135922E-01 1.000000E+21
+ 53 traj.phases.phase.indep_states.states:gam_3 c -1.000000E+21 4.797863E-01 1.000000E+21
+ 54 traj.phases.phase.indep_states.states:gam_4 c -1.000000E+21 4.379029E-01 1.000000E+21
+ 55 traj.phases.phase.indep_states.states:gam_5 c -1.000000E+21 3.873723E-01 1.000000E+21
+ 56 traj.phases.phase.indep_states.states:gam_6 c -1.000000E+21 3.276275E-01 1.000000E+21
+ 57 traj.phases.phase.indep_states.states:gam_7 c -1.000000E+21 2.582142E-01 1.000000E+21
+ 58 traj.phases.phase.indep_states.states:gam_8 c -1.000000E+21 1.789451E-01 1.000000E+21
+ 59 traj.phases.phase.indep_states.states:gam_9 c -1.000000E+21 9.009145E-02 1.000000E+21
+ 60 traj.phases.phase.indep_states.states:gam_10 c -1.000000E+21 -7.423019E-03 1.000000E+21
+ 61 traj.phases.phase.indep_states.states:gam_11 c -1.000000E+21 -1.118970E-01 1.000000E+21
+ 62 traj.phases.phase.indep_states.states:gam_12 c -1.000000E+21 -2.208805E-01 1.000000E+21
+ 63 traj.phases.phase.indep_states.states:gam_13 c -1.000000E+21 -3.314120E-01 1.000000E+21
+ 64 traj.phases.phase.indep_states.states:gam_14 c -1.000000E+21 -4.404248E-01 1.000000E+21
+ 65 traj.phases.phase.indep_states.states:gam_15 c -1.000000E+21 -5.451909E-01 1.000000E+21
+ 66 traj.phases.phase.indep_states.states:gam_16 c -1.000000E+21 -6.436480E-01 1.000000E+21
+ 67 traj.phases.phase.indep_states.states:gam_17 c -1.000000E+21 -7.345244E-01 1.000000E+21
+ 68 traj.phases.phase.indep_states.states:gam_18 c -1.000000E+21 -8.172791E-01 1.000000E+21
+ 69 traj.phases.phase.indep_states.states:gam_19 c -1.000000E+21 -8.919360E-01 1.000000E+21
+ 70 traj.phases.phase.indep_states.states:gam_20 c -1.000000E+21 -9.588926E-01 1.000000E+21
+ 71 traj.phases.phase.indep_states.states:gam_21 c -1.000000E+21 -1.018753E+00 1.000000E+21
+ 72 traj.phases.phase.indep_states.states:gam_22 c -1.000000E+21 -1.072207E+00 1.000000E+21
+ 73 traj.phases.phase.indep_states.states:gam_23 c -1.000000E+21 -1.119950E+00 1.000000E+21
+ 74 traj.phases.phase.indep_states.states:gam_24 c -1.000000E+21 -1.162642E+00 1.000000E+21
+ 75 traj.phases.phase.indep_states.states:v_0 c -1.000000E+21 5.750200E+02 1.000000E+21
+ 76 traj.phases.phase.indep_states.states:v_1 c -1.000000E+21 3.997227E+02 1.000000E+21
+ 77 traj.phases.phase.indep_states.states:v_2 c -1.000000E+21 3.060109E+02 1.000000E+21
+ 78 traj.phases.phase.indep_states.states:v_3 c -1.000000E+21 2.471802E+02 1.000000E+21
+ 79 traj.phases.phase.indep_states.states:v_4 c -1.000000E+21 2.066205E+02 1.000000E+21
+ 80 traj.phases.phase.indep_states.states:v_5 c -1.000000E+21 1.769413E+02 1.000000E+21
+ 81 traj.phases.phase.indep_states.states:v_6 c -1.000000E+21 1.543746E+02 1.000000E+21
+ 82 traj.phases.phase.indep_states.states:v_7 c -1.000000E+21 1.368146E+02 1.000000E+21
+ 83 traj.phases.phase.indep_states.states:v_8 c -1.000000E+21 1.230084E+02 1.000000E+21
+ 84 traj.phases.phase.indep_states.states:v_9 c -1.000000E+21 1.121738E+02 1.000000E+21
+ 85 traj.phases.phase.indep_states.states:v_10 c -1.000000E+21 1.037991E+02 1.000000E+21
+ 86 traj.phases.phase.indep_states.states:v_11 c -1.000000E+21 9.752646E+01 1.000000E+21
+ 87 traj.phases.phase.indep_states.states:v_12 c -1.000000E+21 9.307779E+01 1.000000E+21
+ 88 traj.phases.phase.indep_states.states:v_13 c -1.000000E+21 9.020964E+01 1.000000E+21
+ 89 traj.phases.phase.indep_states.states:v_14 c -1.000000E+21 8.868763E+01 1.000000E+21
+ 90 traj.phases.phase.indep_states.states:v_15 c -1.000000E+21 8.827834E+01 1.000000E+21
+ 91 traj.phases.phase.indep_states.states:v_16 c -1.000000E+21 8.875205E+01 1.000000E+21
+ 92 traj.phases.phase.indep_states.states:v_17 c -1.000000E+21 8.989109E+01 1.000000E+21
+ 93 traj.phases.phase.indep_states.states:v_18 c -1.000000E+21 9.149866E+01 1.000000E+21
+ 94 traj.phases.phase.indep_states.states:v_19 c -1.000000E+21 9.340455E+01 1.000000E+21
+ 95 traj.phases.phase.indep_states.states:v_20 c -1.000000E+21 9.546736E+01 1.000000E+21
+ 96 traj.phases.phase.indep_states.states:v_21 c -1.000000E+21 9.757524E+01 1.000000E+21
+ 97 traj.phases.phase.indep_states.states:v_22 c -1.000000E+21 9.964293E+01 1.000000E+21
+ 98 traj.phases.phase.indep_states.states:v_23 c -1.000000E+21 1.016075E+02 1.000000E+21
+ 99 traj.phases.phase.indep_states.states:v_24 c -1.000000E+21 1.034255E+02 1.000000E+21
+
+ Constraints (i - inequality, e - equality)
+ Index Name Type Lower Value Upper Status Lagrange Multiplier (N/A)
+ 0 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 4.433787E-12 0.000000E+00 9.00000E+100
+ 1 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -3.524292E-12 0.000000E+00 9.00000E+100
+ 2 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -4.433787E-12 0.000000E+00 9.00000E+100
+ 3 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -2.955858E-12 0.000000E+00 9.00000E+100
+ 4 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -3.637979E-12 0.000000E+00 9.00000E+100
+ 5 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -2.273737E-12 0.000000E+00 9.00000E+100
+ 6 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -1.818989E-12 0.000000E+00 9.00000E+100
+ 7 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -1.136868E-12 0.000000E+00 9.00000E+100
+ 8 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -9.094947E-13 0.000000E+00 9.00000E+100
+ 9 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -1.818989E-12 0.000000E+00 9.00000E+100
+ 10 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -4.547474E-13 0.000000E+00 9.00000E+100
+ 11 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 0.000000E+00 0.000000E+00 9.00000E+100
+ 12 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 9.094947E-13 0.000000E+00 9.00000E+100
+ 13 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 9.094947E-13 0.000000E+00 9.00000E+100
+ 14 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 0.000000E+00 0.000000E+00 9.00000E+100
+ 15 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 9.094947E-13 0.000000E+00 9.00000E+100
+ 16 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 -4.547474E-13 0.000000E+00 9.00000E+100
+ 17 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 1.818989E-12 0.000000E+00 9.00000E+100
+ 18 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 4.547474E-13 0.000000E+00 9.00000E+100
+ 19 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 1.818989E-12 0.000000E+00 9.00000E+100
+ 20 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 2.728484E-12 0.000000E+00 9.00000E+100
+ 21 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 2.728484E-12 0.000000E+00 9.00000E+100
+ 22 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 2.728484E-12 0.000000E+00 9.00000E+100
+ 23 traj.phases.phase.continuity_comp.defect_states:r e 0.000000E+00 2.273737E-12 0.000000E+00 9.00000E+100
+ 24 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 -2.205525E-11 0.000000E+00 9.00000E+100
+ 25 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 -1.313083E-11 0.000000E+00 9.00000E+100
+ 26 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 -7.162271E-12 0.000000E+00 9.00000E+100
+ 27 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 -3.637979E-12 0.000000E+00 9.00000E+100
+ 28 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 -1.477929E-12 0.000000E+00 9.00000E+100
+ 29 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 -1.136868E-13 0.000000E+00 9.00000E+100
+ 30 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 0.000000E+00 0.000000E+00 9.00000E+100
+ 31 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 4.547474E-13 0.000000E+00 9.00000E+100
+ 32 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 4.547474E-13 0.000000E+00 9.00000E+100
+ 33 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 4.547474E-13 0.000000E+00 9.00000E+100
+ 34 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 4.547474E-13 0.000000E+00 9.00000E+100
+ 35 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 0.000000E+00 0.000000E+00 9.00000E+100
+ 36 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 3.410605E-13 0.000000E+00 9.00000E+100
+ 37 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 1.136868E-13 0.000000E+00 9.00000E+100
+ 38 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 4.547474E-13 0.000000E+00 9.00000E+100
+ 39 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 5.684342E-13 0.000000E+00 9.00000E+100
+ 40 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 1.477929E-12 0.000000E+00 9.00000E+100
+ 41 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 2.160050E-12 0.000000E+00 9.00000E+100
+ 42 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 2.614797E-12 0.000000E+00 9.00000E+100
+ 43 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 2.785328E-12 0.000000E+00 9.00000E+100
+ 44 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 3.069545E-12 0.000000E+00 9.00000E+100
+ 45 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 3.183231E-12 0.000000E+00 9.00000E+100
+ 46 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 3.467449E-12 0.000000E+00 9.00000E+100
+ 47 traj.phases.phase.continuity_comp.defect_states:h e 0.000000E+00 3.495870E-12 0.000000E+00 9.00000E+100
+ 48 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 1.665335E-15 0.000000E+00 9.00000E+100
+ 49 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 9.992007E-16 0.000000E+00 9.00000E+100
+ 50 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 1.332268E-15 0.000000E+00 9.00000E+100
+ 51 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 1.998401E-15 0.000000E+00 9.00000E+100
+ 52 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 3.164136E-15 0.000000E+00 9.00000E+100
+ 53 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 5.051515E-15 0.000000E+00 9.00000E+100
+ 54 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 6.772360E-15 0.000000E+00 9.00000E+100
+ 55 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 1.018630E-14 0.000000E+00 9.00000E+100
+ 56 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 1.115774E-14 0.000000E+00 9.00000E+100
+ 57 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 1.304339E-14 0.000000E+00 9.00000E+100
+ 58 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 1.312839E-14 0.000000E+00 9.00000E+100
+ 59 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 1.143530E-14 0.000000E+00 9.00000E+100
+ 60 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 8.382184E-15 0.000000E+00 9.00000E+100
+ 61 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 4.940492E-15 0.000000E+00 9.00000E+100
+ 62 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 1.443290E-15 0.000000E+00 9.00000E+100
+ 63 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 -7.771561E-16 0.000000E+00 9.00000E+100
+ 64 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 -2.220446E-15 0.000000E+00 9.00000E+100
+ 65 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 -2.997602E-15 0.000000E+00 9.00000E+100
+ 66 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 -3.108624E-15 0.000000E+00 9.00000E+100
+ 67 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 -3.108624E-15 0.000000E+00 9.00000E+100
+ 68 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 -2.886580E-15 0.000000E+00 9.00000E+100
+ 69 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 -2.442491E-15 0.000000E+00 9.00000E+100
+ 70 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 -2.220446E-15 0.000000E+00 9.00000E+100
+ 71 traj.phases.phase.continuity_comp.defect_states:gam e 0.000000E+00 -1.776357E-15 0.000000E+00 9.00000E+100
+ 72 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 2.518163E-11 0.000000E+00 9.00000E+100
+ 73 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 8.185452E-12 0.000000E+00 9.00000E+100
+ 74 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 3.410605E-12 0.000000E+00 9.00000E+100
+ 75 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 1.705303E-12 0.000000E+00 9.00000E+100
+ 76 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 9.663381E-13 0.000000E+00 9.00000E+100
+ 77 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 6.536993E-13 0.000000E+00 9.00000E+100
+ 78 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 3.694822E-13 0.000000E+00 9.00000E+100
+ 79 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 2.273737E-13 0.000000E+00 9.00000E+100
+ 80 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 4.263256E-14 0.000000E+00 9.00000E+100
+ 81 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 -9.947598E-14 0.000000E+00 9.00000E+100
+ 82 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 -2.273737E-13 0.000000E+00 9.00000E+100
+ 83 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 -3.268497E-13 0.000000E+00 9.00000E+100
+ 84 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 -2.984279E-13 0.000000E+00 9.00000E+100
+ 85 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 -2.984279E-13 0.000000E+00 9.00000E+100
+ 86 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 -2.842171E-13 0.000000E+00 9.00000E+100
+ 87 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 -1.989520E-13 0.000000E+00 9.00000E+100
+ 88 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 -1.136868E-13 0.000000E+00 9.00000E+100
+ 89 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 4.263256E-14 0.000000E+00 9.00000E+100
+ 90 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 8.526513E-14 0.000000E+00 9.00000E+100
+ 91 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 1.705303E-13 0.000000E+00 9.00000E+100
+ 92 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 2.700062E-13 0.000000E+00 9.00000E+100
+ 93 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 3.552714E-13 0.000000E+00 9.00000E+100
+ 94 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 3.694822E-13 0.000000E+00 9.00000E+100
+ 95 traj.phases.phase.continuity_comp.defect_states:v e 0.000000E+00 4.121148E-13 0.000000E+00 9.00000E+100
+ 96 traj.phases.phase->initial_boundary_constraint->ke i 0.000000E+00 4.000000E+00 4.000000E+00 u 9.00000E+100
+
+
+ Exit Status
+ Inform Description
+ 0 Solve Succeeded
+--------------------------------------------------------------------------------
+
False
+
Plotting the results#
+exp_out = traj.simulate()
+
+rad = p.get_val('radius', units='m')[0]
+print(f'optimal radius: {rad} m ')
+mass = p.get_val('size_comp.mass', units='kg')[0]
+print(f'cannonball mass: {mass} kg ')
+angle = p.get_val('traj.phase.timeseries.states:gam', units='deg')[0, 0]
+print(f'launch angle: {angle} deg')
+max_range = p.get_val('traj.phase.timeseries.states:r')[-1, 0]
+print(f'maximum range: {max_range} m')
+
+fig, axes = plt.subplots(nrows=1, ncols=1, figsize=(10, 6))
+
+time_imp = p.get_val('traj.phase.timeseries.time')
+
+time_exp = exp_out.get_val('traj.phase.timeseries.time')
+
+r_imp = p.get_val('traj.phase.timeseries.states:r')
+
+r_exp = exp_out.get_val('traj.phase.timeseries.states:r')
+
+h_imp = p.get_val('traj.phase.timeseries.states:h')
+
+h_exp = exp_out.get_val('traj.phase.timeseries.states:h')
+
+axes.plot(r_imp, h_imp, 'bo')
+
+axes.plot(r_exp, h_exp, 'r--')
+
+axes.set_xlabel('range (m)')
+axes.set_ylabel('altitude (m)')
+
+fig, axes = plt.subplots(nrows=4, ncols=1, figsize=(10, 6))
+states = ['r', 'h', 'v', 'gam']
+for i, state in enumerate(states):
+ x_imp = p.get_val(f'traj.phase.timeseries.states:{state}')
+
+ x_exp = exp_out.get_val(f'traj.phase.timeseries.states:{state}')
+
+ axes[i].set_ylabel(state)
+
+ axes[i].plot(time_imp, x_imp, 'bo')
+ axes[i].plot(time_exp, x_exp, 'r--')
+
+plt.show()
+
/usr/share/miniconda/envs/test/lib/python3.10/site-packages/openmdao/core/group.py:2829: UnusedOptionWarning:'ode_eval' <class ODEEvaluationGroup>: ignored flat_src_indices and/or src_shape because src_indices was not specified.
+/usr/share/miniconda/envs/test/lib/python3.10/site-packages/openmdao/core/group.py:2829: UnusedOptionWarning:'ode_eval' <class ODEEvaluationGroup>: ignored flat_src_indices and/or src_shape because src_indices was not specified.
+/usr/share/miniconda/envs/test/lib/python3.10/site-packages/openmdao/core/group.py:2829: UnusedOptionWarning:'ode_eval' <class ODEEvaluationGroup>: ignored flat_src_indices and/or src_shape because src_indices was not specified.
+
Simulating trajectory traj
+
Done simulating trajectory traj
+optimal radius: 0.04186845473894947 m
+cannonball mass: 2.4194917134095593 kg
+launch angle: 32.01931307429536 deg
+maximum range: 3182.9624372149524 m
+