A Study of Super-Nonlinear Motion of a Simple Pendulum

Haiduke Sarafian

Two identically charged simple pendulums are allowed to swing in a vertical plane, acted upon by gravity and their mutual electrostatic interactive repulsive forces. We show that the inclusion of the electrostatic interaction makes the movement of the pendulums highly nonlinear. To describe the motion, we quantify the relevant time-dependent kinematic and dynamic quantities. Our analysis also includes an extended study of equal and oppositely charged pendulums. Motivated by the outcomes of the calculation, the author manufactured a real-life replica of the study demonstrating the features of the interactive pendulums; a photograph of the replica is included.

Introduction and Motivation

The motion of a simple pendulum under the pull of gravity has been studied for ages. Most standard science and engineering texts have chapters devoted to the analysis of this problem. However, a thorough literature search reveals that the description of the motion of one such pendulum when perturbed by exotic forces is much less extensive. Among many feasible potential scenarios of perturbing the motion of a pendulum in a controlled and quantifiable manner is to cross-pollinate mechanics with electrostatics using a pair of charged pendulums. Based on the results of this analysis, one may argue that the answers to the “what-if” scenarios of this project would have remained unresolved before the advent of the Mathematica era.

This article uses the familiar laws of Newtonian mechanics and electrostatic interactions, as introduced to students in college science and engineering courses. Hence, the physics of this article, its mathematical analysis, and the included Mathematica programs might appeal to a wide range of readers.

To simplify the analysis and to stress the important aspects of the physics of the project, we consider a symmetrical situation, with two identical pendulums.

This is the problem: Consider two identical simple pendulums. Assume each pendulum is composed of a point-mass m, electric charge q, and string length ℓ. Swing the pendulums symmetrically about the vertical reference line that passes through their common pivot and hold them horizontal. Drop them simultaneously and let them swing under gravity. Then analyze the motion of each pendulum in terms of the parameters m, q, and ℓ.

General Strategy

Figure  (left) is a sketch of the problem. It shows the pendulums at an instant where each line makes an angle θ with the vertical reference position. One strategy of analyzing the problem is to identify the forces acting on each mass. Later we consider the alternate energy approach. Because of the symmetry of the assembly, it suffices to consider just one of the point-like masses, say the one on the right. In addition to mechanical forces such as weight (Overscript[w, ->]) and tension (Overscript[T, ->]), we include the electrostatic repulsive force (Subscript[(Overscript[F, ->]), c]). The imposed symmetry confines the orientation of the (Subscript[(Overscript[F, ->]), c]) to the horizontal plane. The analysis in this article can easily be extended to study more challenging asymmetrical cases. Having identified the forces acting on the mass (the free-body diagram), we can apply Newton’s laws of motion. Our strategy is to consider the components of Newton’s second law (Overscript[(Subscript[F, net]), ->])=m (Overscript[a, ->]) along the x and y axes, which gives two separate equations of motion. However, because of the geometrical constraints (i.e., the mass is connected to the pivot with an unstretchable string), the coordinates of the mass x(t) and y(t) are interrelated. This correlation reduces the number of independent equations of motion from two to one. Solving the resulting equation yields the time-dependence of one of the coordinates. The relation between x(t) and y(t) then gives the other coordinate, thus leading to the entire kinematics and dynamics of the problem. The details follow.

Analysis

We begin with Newton’s second law, (Overscript[(Subscript[F, net]), ->])=m (Overscript[a, ->]). Because of the planar motion of the mass, we project this equation along the horizontal and vertical directions. This gives (Overscript[(Subscript[((Subscript[F, net])), x]), ->])=m (Overscript[(Subscript[a, x]), ->]) and (Overscript[(Subscript[((Subscript[F, net])), y]), ->])=m (Overscript[(Subscript[a, y]), ->]). There are advantages to confining the motion of the mass within the first quadrant; we select the downward direction as the positive y axis. These two equations yield

(Subscript[F, c])-T sin θ=m (Overscript[x, ..]),

m g-T cos θ=m (Overscript[y, ..]),

where g is the acceleration of gravity, T is the tension in the line, and (Subscript[F, c])=((k (q^2))/(4 (x^2))) is the electrostatic force, with q being the charge and x the horizontal coordinate of the mass. The accelerations of the mass along the x and y axes are (Overscript[x, ..]) and (Overscript[y, ..]). According to standard notation, (OverDot[x])≡(d/(dt))x, and so on. We rearrange these equations and find their ratio:

((T sin θ)/(T cos θ))=(((Subscript[F, c])-m (Overscript[x, ..]))/(m g-m (Overscript[y, ..]))).

Since, according to Figure , tan θ=(x/y), equation () yields

g-(Overscript[y, ..])=(y/x)(((Subscript[F, c])/m)-(Overscript[x, ..])).

Because the length of the string is constant, we write (x^2)+(y^2)=(ℓ^2), which gives y=+(Sqrt[(ℓ^2)-(x^2)]). Differentiating the latter twice with respect to time yields

-(Overscript[y, ..])=[x (Overscript[x, ..])+((OverDot[x])^2)+(((x (OverDot[x])))^2)((((ℓ^2)-(x^2)))^(-1))]((((ℓ^2)-(x^2)))^(-(1/2))).

Substituting equation () into equation () and simplifying the result gives

(ℓ^2)(x^2)(Overscript[x, ..])+(((x^3)(ℓ^2))/((ℓ^2)-(x^2)))((OverDot[x])^2)+g (x^3)(Sqrt[(ℓ^2)-(x^2)])-((k (q^2))/(4 m))((ℓ^2)-(x^2))=0.

Define a dimensionless variable ξ=(x/ℓ), which gives (OverDot[ξ])=((OverDot[x])/ℓ) and (Overscript[ξ, ..])=((Overscript[x, ..])/ℓ). Substituting these into equation () yields

(ξ^2) (Overscript[ξ, ..])+((ξ^3)/(1-(ξ^2)))((OverDot[ξ])^2)+a (ξ^3)(Sqrt[1-(ξ^2)])-b (1-(ξ^2))=0,

where the two auxiliary constant parameters are a=(g/ℓ) and b=((k (q^2))/(4 m (ℓ^3))). The parameter a is purely mechanical while b is influenced by both mechanics and electrostatics. For pendulums charged with equal and opposite charges and those with no charge, the counter equations of equation () are, respectively,

(ξ^2) (Overscript[ξ, ..])+((ξ^3)/(1-(ξ^2)))((OverDot[ξ])^2)+a (ξ^3)(Sqrt[1-(ξ^2)])+b (1-(ξ^2))=0,

(Overscript[ξ, ..])+(ξ/(1-(ξ^2)))((OverDot[ξ])^2)+a ξ(Sqrt[1-(ξ^2)])=0.

Equations (), (), and () are second-order and highly nonlinear differential equations. For a set of chosen practical values of a and b we apply DSolve; Mathematica fails to produce symbolic solutions. We then apply NDSolve along with the relevant initial conditions, namely, x(t=0)=ℓ and (OverDot[x])(t=0)=0. The solutions are shown in Figure .

values={ℓ->1.0,m->5 (10^(-3)),q-> 2 (10^(-6)),k->9 (10^9),g->9.8};

{a,b}={(g/ℓ),((k (q^2))/(4m (ℓ^3)))}/.values

eqns[n_]:=((ξ[t])^2)ξ''[t]+(((ξ[t])^3)/(1-((ξ[t])^2)))((ξ'[t])^2)+a ((ξ[t])^3) (Sqrt[1-((ξ[t])^2)])+n b (1-((ξ[t])^2))/.values

soleqns=Table[NDSolve[{eqns[n]==0,ξ[1 (10^(-8))]==N[1-(10^(-7)),10],ξ'[1 (10^(-8))]==0},ξ[t],{t,1 (10^(-6)),If[n==-1,3.0,If[n==0,0.6,0.55]]}],{n,-1,1}];

Show[Table[Plot[Evaluate[ξ[t]/.soleqns[[m]]],{t,1 (10^(-6)),If[m==1,3.0,0.6]},AxesLabel->{t,"ξ"},PlotStyle->{Thickness[0.004],GrayLevel[0.15 m],If[m==1,Dashing[{1,0.0001}],Dashing[{0.01,0.001} m]]},PlotRange->{Automatic, {0,1}}],{m,1,3}]]

As one would intuitively expect, in the absence of energy loss (such as friction at the pivot), the electrostatic repulsive force between the two identically charged pendulums would have two distinct impacts. First, it should set the pendulums in steady repetitive oscillations. Second, it prevents the pendulums from touching one another. Both effects are clearly illustrated by the solid curve of Figure . More specifically, the minimum ordinate of the first cycle is the separation distance between the instantly halted pendulums. For oppositely charged pendulums, the attractive electrostatic force causes the pendulums to pull together more strongly than they would otherwise. For the latter as well as the uncharged case, the two pendulums touch one another along the vertical line, ξ=0. The abscissas of the long-dashed and short-dashed curves in Figure  are the fingerprints of these scenarios. The abscissa of the former is less than the latter due to the attraction of the electrostatic force.

Having identified the horizontal coordinate of the pendulum, we evaluate y=+(Sqrt[(ℓ^2)-(x^2)]). Figure  displays the pair {x(t),y(t)}.

Show[Plot[{Evaluate[ξ[t]/.soleqns[[1]]],(Sqrt[1-(((ξ[t]/.soleqns[[1]]))^2)])},{t,1 (10^(-6)),1.2},AxesLabel->{t,Row[{x[t],", ",y[t]}]},PlotStyle->{{{Thickness [0.004],GrayLevel[0]}},{Thickness [0.004],GrayLevel[0.5]}}],GridLines->{Range[0,1.2,0.1],Automatic}]

The curves in Figure  show the features of the coordinates of the pendulum that one would intuitively expect. The black x(t) and gray y(t) curves are completely out of phase. The pendulum begins from its horizontal extreme x(t=0)=1 (the black curve) while its vertical coordinate begins from its zero height extreme y(t=0)=0 (the gray curve). When the pendulum falls, it approaches the vertical reference position. At the instant when it is closest to the other pendulum (the minimum of the black curve), its distance from the vertical position reaches its maximum (the tip of the gray curve).

To confirm that our computation sets the pendulum in a circular path, we apply ParametricPlot to graph its traversed path.

ParametricPlot[Flatten[{Evaluate[ξ[t]/.soleqns[[1]]],Evaluate[-(Sqrt[1-(((ξ[t]/.soleqns[[1]]))^2)])]}],{t,1 (10^(-6)),0.8},AxesLabel->{x[m],y[m]},PlotRange->{{0,1},{0,-1.0}},GridLines->Automatic,PlotStyle->{Thickness[0.004],GrayLevel[0]},ImageSize->250]

Note that {x(t),y(t)} are given numerically, not analytically. Even so, Mathematica lets us evaluate the time derivatives related to kinematic quantities, such as {(OverDot[x])(t),(Overscript[x, ..])(t)} and {(OverDot[y])(t),(Overscript[y, ..])(t)}. We display these in Figure .

{xcoordinate,xspeed,xacc}=Table[D[ξ[t]/.soleqns[[1]],{t,n}],{n,0,2}];
{ycoordinate,yspeed,yacc}=Table[D[(Sqrt[1-(((ξ[t]/.soleqns[[1]]))^2)]),{t,n}],{n,0,2}];

Column[{{
Plot[{xcoordinate,0.5 xspeed,0.03xacc},{t,0.001,3.0},
PlotStyle->{{Thickness[0.004],GrayLevel[0]},
{Thickness [0.004],GrayLevel[0],Dashing[{0.01}]},
{Thickness[0.004],GrayLevel[0],Dashing[{0.02}]}},AxesLabel->{t,Row[{x,", ",(OverDot[x]),", ",(Overscript[x, ".."])}]},PlotRange->All,ImageSize->275],
Plot[{ycoordinate,0.5yspeed,0.08yacc},{t,0.01,1.2},
PlotStyle->{{Thickness[0.004],GrayLevel[0]},
{Thickness [0.004],GrayLevel[0],Dashing[{0.01}]},
{Thickness[0.004],GrayLevel[0],Dashing[{0.02}]}},AxesLabel->{t,Row[{y,", ",(OverDot[y]),", ",(Overscript[y, ".."])}]},PlotRange->All,ImageSize->275]}}]

Graphically speaking, speed is the slope of the position with respect to time; acceleration is the slope of the speed with respect to time. A close inspection of the depicted plots underlines the graphical interrelationships of the corresponding quantities. The author speculates the “noise” in the (Overscript[y, ..])(t) signal (the long-dashed curve) originated from the second-order numerical derivative procedure. The double-hump of (Overscript[y, ..])(t) would have remained undetected had it not been depicted graphically.

An Alternative Approach Using Polar Coordinates

Since each individual mass swings along a circular path, this suggests that the problem might be more efficiently analyzed in polar coordinates.

Earlier, we introduced a dimensionless variable ξ=(x/ℓ). According to Figure , ξ=sin θ. We then evaluate (Overscript[ξ, . ])=(OverDot[θ]) cos θ and (Overscript[ξ, ..])=(Overscript[θ, ..]) cos θ-((OverDot[θ])^2) sin θ. Substituting these into equation () and simplifying yields the equation of motion in the polar coordinate system,

(Overscript[θ, ..])+a sin θ-b ((cos θ)/((sin^2)θ))=0.

Replacing b by -b or 0 results in counterparts of equations () and (), respectively. In the absence of electric charge, b=0 and equation () describes the motion of an uncharged simple pendulum. For small angles, sin θ≃θ and the corresponding linearized equation assumes a textbook form: (Overscript[θ, ..])+(g/ℓ)θ=0. For a charged pendulum b≠0 and, irrespective of the oscillation amplitude, equation () describes the motion of a highly super-nonlinear pendulum. As in the previous section, applying DSolve fails to generate symbolic output, so again we use the same set of initial conditions and solve the equation numerically. The solution of the equation of motion in the polar coordinate system should match the Cartesian results. As shown in Figure , the agreement is perfect.

eqnsθ[n_]:=θ''[t]+a Sin[θ[t]]+n b ((Cos[θ[t]])/((Sin[θ[t]])^2))

soleqnsθ=Table[NDSolve[{eqnsθ[n]==0,θ[1 (10^(-8))]==(π/2),θ'[1 (10^(-8))]==0},θ[t],{t,1 (10^(-8)),3.0}],{n,-1,0}];

Show[Plot[{(180/π)Evaluate[θ[t]/.soleqnsθ[[1]]],(180/π)Evaluate[ArcTan[((ξ[t]/.soleqns[[1]])/(Sqrt[1-(((ξ[t]/.soleqns[[1]]))^2)]))]]},{t,1 (10^(-8)),3.0},AxesLabel->{t,Row[{θ," (deg)"}]},PlotStyle->{{Thickness[0.015],GrayLevel[0.7]},{{Thickness[0.004],GrayLevel[0]}}},PlotRange->{Automatic, {0,90}}]]

It is instructive to visually quantify the impact of the super-nonlinearity of the motion. In Figure , we plot the angular positions of the charged and the uncharged cases. For the uncharged pendulum, for the sake of comparison only, we assume the pendulum oscillates without colliding with the other uncharged pendulum. According to Figure , the impact of nonlinearity is twofold: (1) it decreases the oscillation amplitudes; and (2) it slows the oscillations. The second property is shown by a time-delayed phase shift.

Show[Plot[{(180/π)Evaluate[θ[t]/.soleqnsθ[[2]]],(180/π)Evaluate[ArcTan[((ξ[t]/.soleqns[[1]])/(Sqrt[1-(((ξ[t]/.soleqns[[1]]))^2)]))]]},{t,1 (10^(-8)),3.0},AxesLabel->{t,Row[{θ," (deg)"}]},PlotStyle->{{Thickness[0.01],GrayLevel[0.7]},{{Thickness[0.004],GrayLevel[0]}}},PlotRange->{Automatic, {-90,90}}]]

Knowing the value of θ(t), we evaluate related physical quantities such as the angular speed ω≡(OverDot[θ])(t) and angular acceleration α≡(Overscript[θ, ..])(t). Here again we emphasize the usefulness of Mathematica: it enables us to differentiate a nonanalytic function. Figure  displays the {θ(t),(OverDot[θ])(t),(Overscript[θ, ..])(t)} for b≠0.

repulsive={θcoordinate,θspeed,θacc}=Table[D[θ[t]/.soleqnsθ[[1]],{t,n}],{n,0,2}];

Plot[{(180/π)θcoordinate,(180/π)0.4θspeed,(180/π)0.02θacc},{t,0.001,3.0},
PlotStyle->{{Thickness[0.004],GrayLevel[0]},
{Thickness[0.004],GrayLevel[0],Dashing[{0.01}]},
{Thickness[0.004],GrayLevel[0],Dashing[{0.02}]}},AxesLabel->{t,Row[{θ,", ",(OverDot[θ]),", ",(Overscript[θ, ".."])}]},PlotRange->All]

As we mentioned in the previous section, graphically, the angular speed and the angular acceleration are to be interpreted as the slope of the angular position with respect to time and the slope of the angular speed with respect to time, respectively. The short-dashed and long-dashed curves depicted in Figure  agree with these interpretations.

Phase Diagrams

In the section Analysis, we evaluated a set of kinematic quantities such as θ(t), (OverDot[θ])(t), and (Overscript[θ, ..])(t) to describe the oscillations of the masses. Now, by suppressing the time variable, t, we apply ParametricPlot and graph subsets of these quantities, namely {θ,(OverDot[θ])}, {θ,(Overscript[θ, ..])}, {(OverDot[θ]),(Overscript[θ, ..])}, and {θ,(OverDot[θ]),(Overscript[θ, ..])}. We also display the same sets for the uncharged pendulums. These are all shown in Figure .

neutral={θcoordinate,θspeed,θacc}=Table[D[θ[t]/.soleqnsθ[[2]],{t,n}],{n,0,2}];

s12=Show[{ParametricPlot[Flatten[{neutral[[1]],neutral[[2]]}],{t,1 (10^(-8)),3},PlotStyle->{Thickness[0.004],GrayLevel[0]}],
ParametricPlot[Flatten[{repulsive[[1]],repulsive[[2]]}],{t,1 (10^(-8)),3},PlotStyle->{Thickness[0.004],GrayLevel[0],Dashing[{0.02}]}]},ImageSize->250,AxesLabel->{θ,(OverDot[θ])},AspectRatio->1];

s13=Show[{ParametricPlot[Flatten[{neutral[[1]],neutral[[3]]}],{t,1 (10^(-8)),3},PlotStyle->{Thickness[0.004],GrayLevel[0]}],
ParametricPlot[Flatten[{repulsive[[1]],repulsive[[3]]}],{t,1 (10^(-8)),3},PlotStyle->{Thickness[0.004],GrayLevel[0],Dashing[{0.02}]}]},ImageSize->250,PlotRange->All,AxesLabel->{θ,(Overscript[θ, ".."])},AspectRatio->1];

s23=Show[{ParametricPlot[Flatten[{neutral[[2]],neutral[[3]]}],{t,1 (10^(-8)),3},PlotStyle->{Thickness[0.004],GrayLevel[0]}],
ParametricPlot[Flatten[{repulsive[[2]],0.4repulsive[[3]]}],{t,1 (10^(-8)),3},PlotStyle->{Thickness[0.004],GrayLevel[0],Dashing[{0.02}]}]},ImageSize->250,PlotRange->All,AxesLabel->{(OverDot[θ]),(Overscript[θ, ".."])},AspectRatio->1];

s3D123=Show[{ParametricPlot3D[Flatten[{neutral[[1]],neutral[[2]],0.4neutral[[3]]}],{t,1 (10^(-8)),3},PlotStyle->{Thickness[0.004],GrayLevel[0]}],
ParametricPlot3D[Flatten[{repulsive[[1]],repulsive[[2]],0.4repulsive[[3]]}],{t,1 (10^(-8)),3},PlotStyle->{Thickness[0.004],GrayLevel[0],Dashing[{0.02}]}]},ImageSize->250,PlotRange->All,AxesLabel->{θ,(OverDot[θ]),(Overscript[θ, ".."])},AspectRatio->1,BoxRatios->1];

Show[GraphicsGrid[{{s12,s13},{s23,s3D123}}],ImageSize->350]

With the exception of the top-left graph, a traditional phase diagram, the other three diagrams have seldom, if ever, been discussed in literature. These three graphs are examples demonstrating how Mathematica can be deployed to explore fresh ideas. Their descriptive interpretations are explicit; they are useful graphs assisting our understanding of the physics of the problem. One of the objectives of this article is to demonstrate the impact of the electrostatic interactions of the charged pendulums on the oscillations of the masses. The two scenarios are distinguished from one another by the presence of electric charge. Hence, to avoid the auxiliary potential side effects, such as mechanical collisions of the uncharged masses, we assume the uncharged pendulums pass through each other when they meet. Therefore, the solid curve in the top-left plot of Figure  is symmetrically extended to the θ<0 domain. The charged and the uncharged pendulums begin from the same horizontal position, θ=π/2. Their respective angular position θ and the angular speed (OverDot[θ]) change accordingly. Beyond θ=π/2, the super-nonlinearity of the oscillations causes these two curves to diverge. In addition, the abscissa of the dashed curve is the smallest separation angle of the charged pendulums.

Conclusion

The analysis of the characteristics of perturbed motion of a simple pendulum as presented in this article illustrates the features of nonlinear dynamics and its interface with mechanics and electrostatics. The author’s extensive search of the literature indicated that this analysis is new. The proposed project accomplishes several instructional and research-oriented objectives. The author applied the basic principles of mechanics that are being taught in introductory physics and engineering courses and laid the foundation one step at a time to develop the equations describing the physics of the problem. A quick review of the article will convince the reader that the concept of the project is not hard to grasp; however, its detailed analysis leading to a quantifiable understanding hinges upon the solutions of the challenging equations of motion. Without Mathematica’s powerful and flexible tools, such as NDSolve, as well as its numerical and especially its graphics utilities, the analysis of this project, the super-nonlinear motion of the pendulum, might have remained unnoticed.

Throughout the article, with its analysis and the accompanying Mathematica programs, the author encourages the reader to examine features of extended challenging scenarios. The article thus suggests a road map for explorations of problems in physics and extends the scope of the current status of nonlinear dynamics. For instance, in most texts the nonlinear motion of a simple pendulum is limited to the analysis of “the large amplitude oscillations” of an uncharged, mechanical pendulum. This article extends consideration to additional nonmechanical perturbation forces.

To make the project as comprehensive as possible, the author built a real-life replica of the study. For practical reasons, two 0.4 caliber cylindrical neodymium magnets are used to mimic the effects of the electrostatic repulsive forces; a digital photo of the replica is included in Figure  (right).

Reference

S. Thornton and J. Marion, “Classical Dynamics of Particles and Systems,” 5th ed., Belmont, CA: Brooks/Cole, 2003.