.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples/beginner/plot_betts_10_7.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_examples_beginner_plot_betts_10_7.py: Hypersensitive Control ====================== Objectives ---------- - Show how ``opty`` works with a single differential equation. - Shows how one can improve the accuracy by looking at the solution and changing parameters accordingly. (Here the solution is constant except near the beginning and near the end, so reducing the running time increases the accuracy, without having to resort to more nodes.) - Uses a temporary directory to cache the compiled binary to avoid recompiling if the differential equations do not change. Introduction ------------ This is example 10.7 from [Betts2010]_. As explained there, it was selected to be a very sensitive control problem. **States** - y : state variable **Specifieds** - u : control variable .. GENERATED FROM PYTHON SOURCE LINES 31-37 .. code-block:: Python import numpy as np import sympy as sm import sympy.physics.mechanics as me from opty.direct_collocation import Problem from opty.utils import create_objective_function .. GENERATED FROM PYTHON SOURCE LINES 38-39 Equations of motion. .. GENERATED FROM PYTHON SOURCE LINES 39-46 .. code-block:: Python t = me.dynamicsymbols._t y, u = me.dynamicsymbols('y u') eom = sm.Matrix([-y.diff(t) - y**3 + u]) sm.pprint(eom) .. rst-class:: sphx-glr-script-out .. code-block:: none ⎡ 3 d ⎤ ⎢u(t) - y (t) - ──(y(t))⎥ ⎣ dt ⎦ .. GENERATED FROM PYTHON SOURCE LINES 47-49 The optimization problem and its solution is packed into a function, so that you can easily call it with different parameters. .. GENERATED FROM PYTHON SOURCE LINES 49-86 .. code-block:: Python def solve_optimization(nodes, tf): t0, tf = 0.0, tf num_nodes = nodes interval_value = (tf - t0)/(num_nodes - 1) # Provide some reasonably realistic values for the constants. state_symbols = (y,) specified_symbols = (u,) # Specify the objective function and form the gradient. obj_func = sm.Integral(y**2 + u**2, t) sm.pprint(obj_func) obj, obj_grad = create_objective_function( obj_func, state_symbols, specified_symbols, tuple(), num_nodes, node_time_interval=interval_value, time_symbol=t) # Specify the symbolic instance constraints, as per the example. instance_constraints = ( y.func(t0) - 1, y.func(tf) - 1.5, ) # Create the optimization problem and set any options. prob = Problem(obj, obj_grad, eom, state_symbols, num_nodes, interval_value, instance_constraints=instance_constraints, time_symbol=t, tmp_dir='generated_code') prob.add_option('nlp_scaling_method', 'gradient-based') # Give some rough estimates for the trajectories. initial_guess = np.zeros(prob.num_free) # Find the optimal solution. solution, info = prob.solve(initial_guess) return solution, info, prob .. GENERATED FROM PYTHON SOURCE LINES 87-88 As per the example tf = 10000. .. GENERATED FROM PYTHON SOURCE LINES 88-96 .. code-block:: Python tf = 10000 num_nodes = 501 solution, info, prob = solve_optimization(num_nodes, tf) print(info['status_msg']) print(f'Objective value achieved: {info["obj_val"]:.4f}, as per the book ' f'it is {6.7241}, so the error is: ' f'{(info["obj_val"] - 6.7241)/6.7241*100:.3f} % ') .. rst-class:: sphx-glr-script-out .. code-block:: none ⌠ ⎮ ⎛ 2 2 ⎞ ⎮ ⎝u (t) + y (t)⎠ dt ⌡ b'Algorithm terminated successfully at a locally optimal point, satisfying the convergence tolerances (can be specified by options).' Objective value achieved: 282.5099, as per the book it is 6.7241, so the error is: 4101.453 % .. GENERATED FROM PYTHON SOURCE LINES 97-98 Plot the optimal state and input trajectories. .. GENERATED FROM PYTHON SOURCE LINES 98-100 .. code-block:: Python _ = prob.plot_trajectories(solution) .. image-sg:: /examples/beginner/images/sphx_glr_plot_betts_10_7_001.png :alt: State Trajectories, Input Trajectories :srcset: /examples/beginner/images/sphx_glr_plot_betts_10_7_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 101-102 Plot the constraint violations. .. GENERATED FROM PYTHON SOURCE LINES 102-104 .. code-block:: Python _ = prob.plot_constraint_violations(solution) .. image-sg:: /examples/beginner/images/sphx_glr_plot_betts_10_7_002.png :alt: Constraint violations :srcset: /examples/beginner/images/sphx_glr_plot_betts_10_7_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 105-106 Plot the objective function as a function of optimizer iteration. .. GENERATED FROM PYTHON SOURCE LINES 106-108 .. code-block:: Python _ = prob.plot_objective_value() .. image-sg:: /examples/beginner/images/sphx_glr_plot_betts_10_7_003.png :alt: Objective Value :srcset: /examples/beginner/images/sphx_glr_plot_betts_10_7_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 109-114 With the value of tf = 10000 above, opty converged to a locally optimal point, but the objective value is far from the one given in the book. As per the plot of the solution y(t) it seems, that most of the time y(t) = 0, only at the very beginning and the very end it is different from 0. So, it may make sense to use a smaller tf. Also increasing num_nodes may help. .. GENERATED FROM PYTHON SOURCE LINES 114-122 .. code-block:: Python tf = 8.0 num_nodes = 10001 solution, info, prob = solve_optimization(num_nodes, tf) print(info['status_msg']) print(f'Objective value achieved: {info["obj_val"]:.4f}, as per the book ' f'it is {6.7241}, so the error is: ' f'{(info["obj_val"] - 6.7241)/6.7241*100:.3f} % ') .. rst-class:: sphx-glr-script-out .. code-block:: none ⌠ ⎮ ⎛ 2 2 ⎞ ⎮ ⎝u (t) + y (t)⎠ dt ⌡ b'Algorithm terminated successfully at a locally optimal point, satisfying the convergence tolerances (can be specified by options).' Objective value achieved: 6.7337, as per the book it is 6.7241, so the error is: 0.143 % .. GENERATED FROM PYTHON SOURCE LINES 123-124 Plot the optimal state and input trajectories. .. GENERATED FROM PYTHON SOURCE LINES 124-126 .. code-block:: Python _ = prob.plot_trajectories(solution) .. image-sg:: /examples/beginner/images/sphx_glr_plot_betts_10_7_004.png :alt: State Trajectories, Input Trajectories :srcset: /examples/beginner/images/sphx_glr_plot_betts_10_7_004.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 127-128 Plot the constraint violations. .. GENERATED FROM PYTHON SOURCE LINES 128-130 .. code-block:: Python _ = prob.plot_constraint_violations(solution) .. image-sg:: /examples/beginner/images/sphx_glr_plot_betts_10_7_005.png :alt: Constraint violations :srcset: /examples/beginner/images/sphx_glr_plot_betts_10_7_005.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 131-132 Plot the objective function as a function of optimizer iteration. .. GENERATED FROM PYTHON SOURCE LINES 132-134 .. code-block:: Python _ = prob.plot_objective_value() .. image-sg:: /examples/beginner/images/sphx_glr_plot_betts_10_7_006.png :alt: Objective Value :srcset: /examples/beginner/images/sphx_glr_plot_betts_10_7_006.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 34.757 seconds) .. _sphx_glr_download_examples_beginner_plot_betts_10_7.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_betts_10_7.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_betts_10_7.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_betts_10_7.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_