BIBECHANA Vol. 22, No. 2, August 2025, 116-121 ISSN 2091-0762 (Print), 2382-5340 (Online) Journal homepage: http://nepjol.info/index.php/BIBECHANA Publisher:Dept. of Phys., Mahendra Morang A. M. Campus (Tribhuvan University)Biratnagar From linear to non-linear/chaotic pendulum: a computational study Chhabi Kumar Shrestha1,2, Krishna Raj Adhikari3,∗, Kapil Adhikari4, Narayan Prasad Adhikari2 1Department of Physics, Prithvi Narayan Campus, Tribhuvan University, Pokhara, Nepal 2Central Department of Physics, Tribhuvan University, Kirtipur, Kathmandu, Nepal 3Department of Applied Sciences, Pashchimanchal Campus, Institute of Engineering, Tribhuvan University, Pokhara, Nepal 4Faculty of Science and Technology, Gandaki University, Nepal ∗Corresponding author. Email: adhikarikrishna@wrc.edu.np Abstract In this work, we have used computational techniques to examine how the dynamics of a simple pendulum change from linear to non-linear and chaotic. The graph of phase space, angular displacement versus time, and angular velocity versus time are thoroughly examined in our analysis. A significant shift is observed in these representations, particularly in the graph of angular displacement versus time and angular velocity versus time. As the non-linearity is enhanced, we observe a progressive movement from circular to oval shapes in phase space. In damped and forced pendulum scenarios, similar patterns are observed. In these cases, the graphs display a sinusoidal pattern with a diminishing amplitude with time. Surprisingly, in the phase space of the damped pendulum, a spiral type of graph is observed, demonstrating the intricate relationship between damping effects and non-linearity. This research emphasizes the separatrix’s function as a crucial cutoff point where the pendulum’s motion changes from linear to chaotic. Keywords Phase space; Separatrix; Damped pendulum; Trajectory; Runge-Kutta method; Harmonic oscillator Article information Manuscript received: December 12, 2024; Revised: January 7, 2025; Accepted: January 10, 2025 DOI https://doi.org/10.3126/bibechana.v22i2.72570 This work is licensed under the Creative Commons CC BY-NC License. https://creativecommons. org/licenses/by-nc/4.0/ 1 Introduction Galileo was the first to observe that the swinging time of a lamp in a cathedral remained consistent regardless of its swing’s size, at least for the small swings he could observe [1]. Huygens invented the first pendulum clock in 1656 [2], which was a big step forward in timekeeping from earlier techniques. 116 http://nepjol.info/index.php/BIBECHANA adhikarikrishna@wrc.edu.np https://doi.org/10.3126/bibechana.v22i2.72570 https://creativecommons.org/licenses/by-nc/4.0/ https://creativecommons.org/licenses/by-nc/4.0/ Chhabi Kumar Shrestha et al./ BIBECHANA 22 (2025) 116-121 117 The pendulum was thus the first oscillator of signif- icant practical significance. The study of pendulum dynamics plays an important role in classical me- chanics, and simple harmonic motion (SHM) is well understood for small oscillations. If the particle in the case of periodic motion moves back and forth, then the motion is called oscillatory or vibratory motion. The harmonic motion of the simplest type, having constant amplitude and single frequency, is called simple harmonic motion. The body is said to be in simple harmonic motion (SHM) if the ac- celeration of the body is directly proportional to the displacement and is directed always toward the mean position. If ’a’ be the acceleration of the body and ’x’ be the displacement of the body from the mean position, then acceleration ∝ displacement i.e., a ∝ x ⇒ a = −kx where ’k’ is a constant, called a force constant or spring constant, and negative sign indicates accel- eration is always directed opposite to the displace- ment or motion of an object [3]. Suppose a particle of mass ’m’ executing SHM. If ’x’ be the displacement of the particle from the equilib- rium position at any instant of time ’t’, then from Hooke’s law Restoringforce(F ) ∝ x ⇒ F = −kx ⇒ d2x dt2 = −kx m = −w2x Where w2 = k m The general solution of this equa- tion is x = A1e iw0t +A2e −iw0t Phase space is the momentum, and position coor- dinates are set to define the dynamic system. The idea that unifies quantum and classical mechanics is crucial to understanding physics. The space of all possible states for a physical system is known as phase space in classical mechanics. Friedrich Henri Poincare, Ludwig Boltzmann, and Josiah Willard Gibbs invented the idea of phase space in the late 1800 s [4]. A dynamical system’s state varies with time, as represented by a point in phase space, and this is known as the phase trajectory. A phase plane that shows the angle θ and the angular ve- locity dθ dt provides a clear visual representation of the dynamics of a pendulum. Certain curves or points, referred to as “attractors", are important in this graphical representation. No matter where the motion begins, attractors are crucial because they symbolize the final states that all potential pendu- lum trajectories can converge toward over time. In the case of real-world systems, when large angu- lar displacement and the damping term are taken into account, a simple pendulum exhibits non-linear and chaotic behavior. While many studies have an- alyzed these behaviors separately, a comprehensive computational study that systematically examines the transition from linear motion to non-linearity and chaos is still lacking. This work has employed the numerical technique Runge-Kutta method to investigate chaotic behavior in pendulums, but fur- ther exploration is required to fill the gap between classical and modern computational approaches. This study aims to address this gap by employ- ing a computational framework to analyze phase- space trajectories and quantify the transition from periodic to chaotic motion. Understanding these transitions has important applications in fields such as seismology, robotics, and precision timekeeping, where chaotic oscillations impact system stability and performance. 2 Materials and Methods A simple pendulum is a heavy-point object sus- pended by a flexible, weightless, and inextensible string that can oscillate freely in a vertical plane. When a bob of a simple pendulum is displaced through a small angle from its mean position, then it begins to vibrate around its mean position. Such motion of the point mass is called simple harmonic motion. Suppose a bob of a simple pendulum of mass ‘m’ whose effective length is ‘l’. Let the bob be displaced through ‘x’ such that the angle θ is very small. Then various forces acting on the bob are: 1. weight ‘mg’ acting vertically downward. 2. Tension ’T’ acting on a string towards the point of suspension Here, weight ‘mg’ can be resolved into two compo- nents. The component mgcosθ balances tension T, and another component mgsinθ is restoring force, which causes the oscillation of the simple pendu- lum. According to Newton’s second law, F = ma = m d2x dt2 m d2x dt2 = −mgsinθ The negative sign indicates that restoring force and displacement are in opposite directions. If we have x = lθ then above equation becomes ml d2θ dt2 = −mgsinθ d2θ dt2 = −g L sinθ (1) Chhabi Kumar Shrestha et al./ BIBECHANA 22 (2025) 116-121 118 Using Maclaurin’s theorem, we get, d2θ dt2 = −g L (θ − θ3 3! + θ5 5! − ......) (2) Taking the first term of the right-hand side, we get, θ̈ = −g L θ (3) Which has the solution of the form θ(t) = θ0cos √ g L t (4) Its time period is T = 2π √ L g [6] 2.1 Damped harmonic oscillator A harmonic oscillator is called a damped harmonic oscillator in which the oscillations are damped on account of a resistive or damping force, with its am- plitude progressively decreasing to zero. When a simple harmonic oscillator is subjected to a damp- ing force that is proportional to its velocity and an external periodic force, then we write an equation of its motion [7]. d2x dt2 + γ dx dt + ω2 0x = a0sinwt (5) Where F (t) = ma0sinwt is the driving force with angular frequency ω.In the absence of driving force i.e., at a = 0, the real solutions of equation (7)at t → ∞ are: x0(t) = c1e −(γ+ √ γ2−4ω2 0)t/2+c2e −(γ− √ γ2−4ω2 0)t/2, γ2 − 4w2 0 > 0 (6) x0(t) = c1e −γt/2 + c2e −γt/2, γ2 − 4w2 0 = 0 (7) x0(t) = c1e −γt/2cos (√ −γ2 + 4ω2 0t/2 ) + c2e −γt/2cos (√ −γ2 + 4ω2 0t/2 ) γ2 − 4w2 0 < 0 (8) If a0 > 0, then the general solution is obtained from the sum of the special solution, xs(t), and the solu- tion of homogeneous equations is x0(t). The special solution is obtained as xs(t) = Asinwt+Bcoswt1001[8] Putting this equation in equation (7), we get, xx(t) = a0[(w 2 0 − w2)coswt+ γwsinwt] (w2 0 − w2)2 + γ2w2 (9) So, the general solution is given by x(t) = x0(t) + xs(t) (10) In the case of resonance without damping w = w0 and γ = 0. In this case, the solution becomes x(t) = c1coswt+c2sinwt+ a0 4w2 (coswt+2wtsinwt) (11) 2.2 Forced damped pendulum Let us consider a simple pendulum in a uniform gravitational field in which its motion is damped by the force that is proportional to its velocity and is under the action of the vertical harmonic exter- nal driving force then its equation of motion be- comes [9]. d2θ dt2 + γ dθ dt + w2 0sinθ = −2Acoswtsinθ (12) Where θ is the angle made by the pendulum with the vertical, γ is the damping coefficient, w2 0 = g/l is the natural angular frequency of the simple pen- dulum, ω is the angular frequency of driving force. The 2A is the amplitude of the oscillation. 3 Methodology The current study is both computational and the- oretical in nature. It mainly focuses on the phase- space trajectory of the chaotic pendulum. Using the Runge-Kutta method, a numerical solution is obtained. Fortran 90 is used to run the file, and Gnuplot is used to plot the phase space. Fortran is used to program scientific and mathematical appli- cations. The name FORTRAN is an acronym for Formula Translation. It is a high-level program- ming language. There have been different versions of Fortran since the 1950s when work on it began at IBM. The equation (1) can be solved by using the Runge-Kutta (RK) method [10, 11].From the pedagogical perspective, Fortran90 offers a basic method for numerical computing. It is a well-known language for high-performance numerical compu- tation that is very effective at solving differen- tial equations involving large data sets. Fortran is based on a numerical algorithm that aids re- searchers and students in expanding their under- standing of numerical analysis, whereas Python and MATLAB are based on built-in libraries. Fortran 90 is an excellent tool for solving complicated phys- ical models, such as chaotic pendulums, because it offers a better-optimized solution for iterative nu- merical computations. The simulation parameters Chhabi Kumar Shrestha et al./ BIBECHANA 22 (2025) 116-121 119 employed in this investigation are the pendulum’s length of 2.0 meters, damping coefficient of 0.5, sim- ulation time step of 30 seconds, and angular velocity of 3.1 rad/sec. Even though there are many numer- ical techniques, such as Euler’s method and Verlet integration, we have chosen to employ the Runge- Kutta method (RK) because it offers a higher-order accurate solution and is frequently applied to non- linear dynamical systems because of its enhanced stability. The first-order Euler’s approach makes it less precise, requiring a very small time step to minimize errors. It is not dependable for long-term simulations due to its significant truncation errors. Likewise, Verlet integration works well in molecu- lar dynamics simulations but is less appropriate for chaotic differential equations. 4 Results and Discussions In this study entitled “From Linear to Non- linear/Chaotic Pendulum: A Computational Study," we investigated the changes in a pendu- lum system’s dynamic behavior. Considering the equation θ̈ = −g L sinθ in the beginning, we tracked the pendulum’s angular displacement θ against time, which is shown in Figure 2 for the angle of oscillation 90.In the linear regime—where the an- gular displacement remains small—the pendulum displayed simple harmonic motion, which is char- acterized by smooth, sinusoidal oscillations. As the amplitude increases, i.e., the angle of oscillation is increased and moved away from basic harmonic behavior, the motion becomes more complex and enters into the non-linear zone. Chaotic dynam- ics began as the angle grew closer to and beyond 180, at which point the non-linear behavior became more noticeable. We then looked at the relationship between an- gular velocity and time, which is shown in Figure 3. A constant amplitude was maintained by the harmonic oscillation of the angular velocity in the linear regime. Significant amplitude and frequency changes were observed as the system entered the non-linear region. Again considering the equation d2θ dt2 +γ dθ dt +w2 0sinθ = −2Acoswtsinθ, we have plot- ted the graph of θ versus time, θ̇ versus time and θ̇ versus θ. When the damping term was present, additional complexity was created because energy dissipation changed the oscillatory patterns, which, over time, caused the motion to get slower grad- ually. Plots showing the relationship between θ̈ versus θ the phase space analysis, shown in Figure 4, helps us to understand the behavior of the sys- tem over time. For small oscillations, the phase-space trajec- tory forms closed loops, confirming agreement with theoretical expectations of simple harmonic motion. As the amplitude increases, the deviation from the simple harmonic motion is significant, leading to distortion in the phase space indicating the chaotic behavior. The separatrix is the clear indicator for such transition, marking the boundary between reg- ular and chaotic motion and observed in the phase space plots beyond 1800. These findings are also in agreement with the theoretical models. When the pendulum enters the chaotic regime, the phase- space diagram exhibits irregular, scattered trajecto- ries, which is consistent with known chaotic behav- ior in dynamical systems. Similar trends were ob- served in the case of the damped pendulum, where damping affected the rate of energy dissipation in the system shown in Figure 5, and Figure 6. Plots of the damped pendulum’s phase space revealed spi- ral paths that eventually converged to fixed points, signifying the end of the motion of the system as depicted in Figure 7. Figure 2: Graph of θ versus t for θ = 90 and dt=0.0 Similarly, the graph of angular velocity θ̇ versuss time is shown in Figure 3. Figure 3: Graph of θ̇ versus t for θ = 90 and dt=0.0 The phase space of the nonlinear pendulum is Chhabi Kumar Shrestha et al./ BIBECHANA 22 (2025) 116-121 120 shown in Figure 4. Figure 4: Graph of θ̇ versus t for dt =0.02sec Figure 5: Graph of θ versus t for γ = 0.5 , dt=0.01,ω = 3.1, ω0 = 3.1 Figure 6: Graph of θ̇ versus t for γ = 0.5 , dt=0.01,ω = 3.1, ω0 = 3.1 for γ = 0.5 , dt=0.01, ω = 3.1, ω0 = 3.1 Figure 7: Graph of θ̇ versus θ 5 Conclusions This study offers insightful information about how linear to non-linear behavior changes in pendulum systems. It clarifies the complex dynamics un- derlying pendulum motion and advances our un- derstanding of chaotic behavior by fusing rigor- ous analysis with numerical simulations. The find- ings have implications beyond theoretical issues and may be used in engineering, physics, and other dis- ciplines. The study of non-linear dynamics and its applications to complex systems requires further in- vestigation. 6 Acknowledgments We extend my sincere gratitude to Saddam Hussein Dhobi and Sandip Baral for their invaluable guid- ance, unwavering support, and insightful feedback throughout the entire duration of this work, “From Linear to Non-linear/Chaotic Pendulum: A Com- putational Study." One of the authors, Chhabi Ku- mar Shrestha, acknowledges the UGC, Nepal, for the PhD grant award: PhD-80/81-S&T-09. References [1] Laura Fermi and Gilberto Bernardini. Galileo and the scientific revolution. Courier Corpora- tion, 2003. [2] Filip Buyse. Galileo, huygens and the pen- dulum clock: Isochronism and synchronicity1. Societate si politica, 11(2):5–12, 2017. [3] D. S. Mathur. Mechanics. S. Chand Publish- ing, 2007. [4] David D Nolte. The tangled tale of phase space. Physics today, 63(4):33–38, 2010. Chhabi Kumar Shrestha et al./ BIBECHANA 22 (2025) 116-121 121 [5] Neha Aggarwal, Nitin Verma, and P Arun. Simple pendulum revisited. European journal of physics, 26(3):517, 2005. [6] P. L. DeVries and R. P. Wolf. A first course in computational physics. Computers in Physics, 8(2):178–179, 1994. [7] Konstant. Computational physics. http://www.physics.ntua.gr/~konstant/ ComputationalPhysics. [8] K. N. Anagnostopoulos. Computational Physics, Vol I: A Practical Introduction to Computational Physics and Scientific Comput- ing. Konstantinos Anagnostopoulos, 2014. [9] M. Suzuki and I. S. Suzuki. Physics of simple pendulum. 2019. [10] D. W. Zingg and T. T. Chisholm. Runge–kutta methods for linear ordinary differential equa- tions. Applied Numerical Mathematics, 31(2):227–238, 1999. [11] A. Tocino and R. Ardanuy. Runge–kutta methods for numerical solution of stochastic differential equations. Journal of Computa- tional and Applied Mathematics, 138(2):219– 241, 2002. http://www.physics.ntua.gr/~konstant/ComputationalPhysics http://www.physics.ntua.gr/~konstant/ComputationalPhysics Introduction Materials and Methods Damped harmonic oscillator Forced damped pendulum Methodology Results and Discussions Conclusions Acknowledgments