.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "api/gallery/classical/newtonian/plot_01_cannonball.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_api_gallery_classical_newtonian_plot_01_cannonball.py: Free fall with an oblique velocity: the cannonball problem ============================================================= :class:`~physicskit.classical.systems.newtonian.ProjectileMotion` is free fall under constant gravity :math:`g` with an oblique initial velocity -- the historical "cannonball" problem. Without drag it is the separable Hamiltonian .. math:: H(q, p) = \frac{|p|^2}{2m} + m g y, \qquad q = (x, y) , whose exact solution is Galileo's parabolic trajectory :math:`x(t) = v_{x0} t,\ y(t) = y_0 + v_{y0} t - \tfrac12 g t^2`. With ``drag_coeff`` :math:`> 0` a quadratic air-resistance force :math:`F_\text{drag} = -c\,|v|\,v` is added, making the system genuinely dissipative (energy decreases monotonically; there is no closed-form trajectory). Fires a projectile at 45 degrees and shows: 1. Its trajectory next to Galileo's closed-form parabola -- since the force is constant, the symplectic Yoshida4 integrator reproduces the exact solution to floating-point precision, at any step size. 2. The effect of quadratic air drag: shorter range, steeper descent, and monotonically decreasing energy (drag is genuinely non-conservative, so it needs ``method="rk4"`` -- the symplectic backends only see the conservative gravity force and would silently ignore ``drag_coeff``). .. GENERATED FROM PYTHON SOURCE LINES 32-37 .. code-block:: Python import matplotlib.pyplot as plt import numpy as np from physicskit.classical.systems.newtonian import ProjectileMotion .. GENERATED FROM PYTHON SOURCE LINES 38-43 Vacuum trajectory vs. the analytic parabola --------------------------------------------- Because the force is constant, a symplectic integrator reproduces Galileo's exact closed-form parabola to floating-point precision, at any step size -- not merely "small energy drift." .. GENERATED FROM PYTHON SOURCE LINES 43-68 .. code-block:: Python proj = ProjectileMotion.from_launch(speed=20.0, angle_deg=45.0, g=9.81) result = proj.integrate((0, 2.9), dt=1e-3, method="yoshida4") x_an, y_an = ProjectileMotion.analytic_trajectory(20.0, 45.0, 9.81, result.t) rng, h_max = ProjectileMotion.range_and_max_height(speed=20.0, angle_deg=45.0, g=9.81) print(f"closed-form range = {rng:.3f} m, max height = {h_max:.3f} m") print(f"max |numeric - analytic| position error: {np.max(np.abs(result.q[:, 1] - y_an)):.2e} m") fig1, axes = plt.subplots(1, 2, figsize=(9, 4)) axes[0].plot(result.q[:, 0], result.q[:, 1], color="steelblue", lw=2, label="numeric (Yoshida4)") axes[0].plot(x_an, y_an, "--", color="orange", lw=1.2, label="analytic") axes[0].axhline(0, color="0.7", lw=0.8) axes[0].set_xlabel("x (m)") axes[0].set_ylabel("y (m)") axes[0].set_title("Trajectory") axes[0].legend() axes[0].set_aspect("equal") axes[1].plot(result.t, np.abs(result.q[:, 1] - y_an), color="crimson") axes[1].set_xlabel("t (s)") axes[1].set_ylabel("|y_numeric - y_analytic| (m)") axes[1].set_title("Numeric vs. closed form") fig1.tight_layout() .. image-sg:: /api/gallery/classical/newtonian/images/sphx_glr_plot_01_cannonball_001.png :alt: Trajectory, Numeric vs. closed form :srcset: /api/gallery/classical/newtonian/images/sphx_glr_plot_01_cannonball_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none closed-form range = 40.775 m, max height = 10.194 m max |numeric - analytic| position error: 3.33e-12 m .. GENERATED FROM PYTHON SOURCE LINES 69-79 Vacuum vs. quadratic air drag --------------------------------- Real cannonballs feel quadratic air drag, :math:`F_\text{drag} = -k |v| v`. This depends on velocity (momentum), so it cannot enter a separable ``force(q, t)`` callback -- it is invisible to Verlet/Yoshida4, which is why :class:`~physicskit.classical.systems.newtonian.ProjectileMotion` *refuses* ``method="yoshida4"``/``"verlet"`` once ``drag_coeff`` is nonzero, rather than silently dropping the drag term. Use ``method="rk4"`` instead, where energy now decreases monotonically -- the whole point of drag. .. GENERATED FROM PYTHON SOURCE LINES 79-98 .. code-block:: Python vac = ProjectileMotion.from_launch(speed=20.0, angle_deg=45.0, g=9.81) drag = ProjectileMotion.from_launch(speed=20.0, angle_deg=45.0, g=9.81, drag_coeff=0.05) res_vac = vac.integrate((0, 2.9), dt=1e-3, method="yoshida4") res_drag = drag.integrate((0, 2.9), dt=1e-3, method="rk4") print(f"vacuum final energy: {res_vac.energy[-1]:.3f} J (started at {res_vac.energy[0]:.3f} J)") print(f"drag final energy: {res_drag.energy[-1]:.3f} J (started at {res_drag.energy[0]:.3f} J, monotonic decrease)") fig2, ax = plt.subplots(figsize=(6, 5)) ax.plot(res_vac.q[:, 0], res_vac.q[:, 1], color="steelblue", lw=2, label="vacuum") ax.plot(res_drag.q[:, 0], res_drag.q[:, 1], color="firebrick", lw=2, label="quadratic drag") ax.axhline(0, color="0.7", lw=0.8) ax.set_xlabel("x (m)") ax.set_ylabel("y (m)") ax.set_title("Drag shortens the range and steepens the descent") ax.legend() ax.set_aspect("equal") .. image-sg:: /api/gallery/classical/newtonian/images/sphx_glr_plot_01_cannonball_002.png :alt: Drag shortens the range and steepens the descent :srcset: /api/gallery/classical/newtonian/images/sphx_glr_plot_01_cannonball_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none vacuum final energy: 200.000 J (started at 200.000 J) drag final energy: 1.723 J (started at 200.000 J, monotonic decrease) .. GENERATED FROM PYTHON SOURCE LINES 99-102 See also :class:`~physicskit.classical.systems.newtonian.KeplerSystem`, which is the same "central conservative force, integrated symplectically" idea, generalized from uniform gravity down to an inverse-square field. .. GENERATED FROM PYTHON SOURCE LINES 102-104 .. code-block:: Python plt.show() .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.819 seconds) .. _sphx_glr_download_api_gallery_classical_newtonian_plot_01_cannonball.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_01_cannonball.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_01_cannonball.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_01_cannonball.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_