phase portrait - Nonlinear damped pendulum

[Giovanni_62]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42

import numpy as np import matplotlib.pyplot as plt from scipy.integrate import solve_ivp   # Parameters g = 9.81       # Acceleration due to gravity (m/s^2) L = 1.0        # Length of the pendulum (m) b = 0.1        # Damping coefficient   # Define the ODEs def pendulum_ode(t, state):     theta, omega = state     dtheta_dt = omega     domega_dt = - (g / L) * np.sin(theta) - b * omega     return [dtheta_dt, domega_dt]   # Create a grid of initial conditions in the phase space theta_range = np.linspace(-2 * np.pi, 2 * np.pi, 3) omega_range = np.linspace(-10, 10, 7)   # Initialize lists to store trajectories trajectories = []   # Integrate the ODEs for a subset of initial conditions for theta_0 in theta_range:     for omega_0 in omega_range:         initial_state = [theta_0, omega_0]         t_span = (0, 10)  # Time span for integration         sol = solve_ivp(pendulum_ode, t_span, initial_state, t_eval=np.linspace(*t_span, 1000))         trajectories.append(sol.y)   # Plot the phase portrait with restricted x and y limits plt.figure(figsize=(8, 6)) for traj in trajectories:     plt.plot(traj[0], traj[1], linewidth=0.5) plt.xlim(-2 * np.pi, 2 * np.pi)  # Restrict x-axis limits plt.ylim(-10, 10)  # Restrict y-axis limits plt.xlabel('Angle (radians)') plt.ylabel('Angular Velocity (radians/s)') plt.title('Phase Portrait of Damped Nonlinear Pendulum') plt.grid(True) plt.show()

Someone help me improve the program, so that we have a visualization of the various trajectories, for example as in

Plot phase diagram

[deanhystad]

Looks interesting:

https://stackoverflow.com/questions/5870...oordinates

http://be150.caltech.edu/2018/handouts/l...raits.html

https://pysd-cookbook.readthedocs.io/en/...raits.html