.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples\plot_betts_10_148_opty.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_plot_betts_10_148_opty.py: Zermelo's Problem (betts-10-148) ================================ This is example 10.148 from *John T. Betts, Practical Methods for Optimal Control Using NonlinearProgramming*, 3rd edition, Chapter 10: Test Problems. The goal is to minimize the final time :math:`t_f` to reach the point (0 / 0), having started at the point (3.5 / -1.8). More information about Zermelo's navigation problem can be found at https://en.wikipedia.org/wiki/Zermelo%27s_navigation_problem **States** - :math:`x, y` : state variables **Controls** - :math:`\theta` : control variable .. GENERATED FROM PYTHON SOURCE LINES 23-31 .. code-block:: Python import numpy as np import sympy as sm import sympy.physics.mechanics as me import time from opty import Problem from opty.utils import MathJaxRepr .. GENERATED FROM PYTHON SOURCE LINES 32-33 Equations of motion. .. GENERATED FROM PYTHON SOURCE LINES 33-47 .. code-block:: Python t = me.dynamicsymbols._t x, y = me.dynamicsymbols('x, y') theta = me.dynamicsymbols('theta') h = sm.symbols('h') V, c = sm.symbols('V, c') eom = sm.Matrix([ -x.diff(t) + V*sm.cos(theta) + c*y, -y.diff(t) + V*sm.sin(theta), ]) MathJaxRepr(eom) .. raw:: html
$$\left[\begin{matrix}V \cos{\left(\theta \right)} + c y - \dot{x}\\V \sin{\left(\theta \right)} - \dot{y}\end{matrix}\right]$$


.. GENERATED FROM PYTHON SOURCE LINES 48-50 Define and Solve the Optimization Problem ----------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 50-60 .. code-block:: Python num_nodes = 2001 t0 = 0*h tf = (num_nodes - 1) * h interval_value = h state_symbols = (x, y) unkonwn_input_trajectories = (theta, ) par_map = {V: 1.0, c: -1.0} .. GENERATED FROM PYTHON SOURCE LINES 61-62 Specify the objective function and form the gradient. .. GENERATED FROM PYTHON SOURCE LINES 62-74 .. code-block:: Python def obj(free): return free[-1] def obj_grad(free): grad = np.zeros_like(free) grad[-1] = 1.0 return grad .. GENERATED FROM PYTHON SOURCE LINES 75-77 Set the instance constraints. Forcing :math:`h \geq 0` will avoid physically meaninigless negative time intervals. .. GENERATED FROM PYTHON SOURCE LINES 77-88 .. code-block:: Python instance_constraints = ( x.func(t0) - 3.5, y.func(t0) + 1.8, x.func(tf), y.func(tf), ) bounds = { h: (0.0, 0.5), } .. GENERATED FROM PYTHON SOURCE LINES 89-90 Create the optimization problem and set any options. .. GENERATED FROM PYTHON SOURCE LINES 90-103 .. code-block:: Python prob = Problem( obj, obj_grad, eom, state_symbols, num_nodes, interval_value, instance_constraints=instance_constraints, known_parameter_map=par_map, bounds=bounds, time_symbol=t, ) .. GENERATED FROM PYTHON SOURCE LINES 104-105 Give some rough estimates for the trajectories. .. GENERATED FROM PYTHON SOURCE LINES 105-107 .. code-block:: Python initial_guess = np.ones(prob.num_free) * 0.1 .. GENERATED FROM PYTHON SOURCE LINES 108-109 Find the optimal solution. .. GENERATED FROM PYTHON SOURCE LINES 109-119 .. code-block:: Python start_time = time.time() solution, info = prob.solve(initial_guess) end_time = time.time() print(f"Solving time: {end_time - start_time:.2f} seconds") print(info['status_msg']) Jstar = 5.26493205 print(f"Objective value achieved: {(num_nodes-1)*info['obj_val']:.4f} ", f"as per the book it is {Jstar:.4f}, so the error is: " f"{((num_nodes-1)*info['obj_val'] - Jstar)/Jstar*100:.3f} % ") .. rst-class:: sphx-glr-script-out .. code-block:: none Solving time: 3.42 seconds b'Algorithm terminated successfully at a locally optimal point, satisfying the convergence tolerances (can be specified by options).' Objective value achieved: 5.2637 as per the book it is 5.2649, so the error is: -0.023 % .. GENERATED FROM PYTHON SOURCE LINES 120-121 Plot the optimal state and input trajectories. .. GENERATED FROM PYTHON SOURCE LINES 121-123 .. code-block:: Python _ = prob.plot_trajectories(solution) .. image-sg:: /examples/images/sphx_glr_plot_betts_10_148_opty_001.png :alt: State Trajectories, Input Trajectories :srcset: /examples/images/sphx_glr_plot_betts_10_148_opty_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 124-125 Plot the constraint violations. .. GENERATED FROM PYTHON SOURCE LINES 125-127 .. code-block:: Python _ = prob.plot_constraint_violations(solution) .. image-sg:: /examples/images/sphx_glr_plot_betts_10_148_opty_002.png :alt: Constraint violations :srcset: /examples/images/sphx_glr_plot_betts_10_148_opty_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 128-129 Plot the objective function as a function of optimizer iteration. .. GENERATED FROM PYTHON SOURCE LINES 129-130 .. code-block:: Python _ = prob.plot_objective_value() .. image-sg:: /examples/images/sphx_glr_plot_betts_10_148_opty_003.png :alt: Objective Value :srcset: /examples/images/sphx_glr_plot_betts_10_148_opty_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 11.363 seconds) .. _sphx_glr_download_examples_plot_betts_10_148_opty.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_148_opty.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_betts_10_148_opty.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_betts_10_148_opty.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_