3. Oscillating Systems

Legal Notice, Disclaimer, and Book Ordering Information:

Null

Setup

\(Clear["Global`*`"]; 
Off[General::spell]; 
Off[General::spell1]; 
    \)

3.3 Problems

3.3.2  Nonlinear Oscillations

Problem 1:  Duffing's Oscillator Equation

<< "D:\Mma\Inits\ch3init.m"

Clear["Global`*`"];

Part (a)   Illustrate the double-well and single-well behavior of the potential by plottingthe potential as a function of a and q.  Pick values of a that correspond to the two kinds of behavior and plot the potential as a function of q.  Plot the phase diagram for the double-well and single-well behavior.

Fke = q''[t]

Fpe = -a q[t]

Fpert  =  b q[t]^3

force = Fke + Fpe + Fpert

T = Integrate[q'[t] Fke, t]

V0 = Integrate[q'[t] Fpe, t]

V1 = Integrate[q'[t] Fpert, t]

energy = T + V0 + V1

\(\(Plot3D[V0 + V1 //. {q[t] -> q, b -> 0.05} // Evaluate

      	, {q, \(-15\), 15}
\ \ , {a, \(-1\), 4}
\ \ , 
      AxesLabel -> {"q", "a", "Potential"}
	, 
      BoxRatios -> {1, 1, 1}
	, PlotRange -> {\(-100\), 100}
	, 
      ViewPoint -> {0, 3, 2}]; \)\)

\(\(plot1 = 
      Plot[V0 + V1 //. {q[t] -> q, b -> 0.05, a -> \(+4\)}

        	, {q, \(-15\), 15}
	, PlotLabel -> "Two-Well Potential"]; 
    \)\)

values = {q[t] -> q, q'[t] -> v, b -> 0.05}

\(\(cplot1 = 
      ContourPlot[\(energy //. values\) //. {a -> \(-1\)} // Evaluate

        	, {q, \(-15\), 15}
	, {v, \(-15\), 15}
	, Contours -> 8
	, 
        PlotPoints -> 20
	, ColorFunction -> Hue]; \)\)

\(\(cplot2 = 
      ContourPlot[\(energy //. values\) //. {a -> \(+4\)} // Evaluate

        	, {q, \(-15\), 15}
	, {v, \(-15\), 15}
	, Contours -> 8
	, 
        PlotPoints -> 20
	, ColorFunction -> Hue]; \)\)

Show[GraphicsArray[{cplot1, cplot2}], 
      DisplayFunction -> $DisplayFunction];

Part (b)  Numerically solve for the trajectories and plot q[t] and q'[t].

\(eq1 = {\((force /. {\(\(q'\)'\)[t] -> \(p'\)[t]})\) == 0
, 
        \(q'\)[t] == p[t]}\)

initial = {p[0] == 0.001, q[0] == 0}

eq2 = Join[eq1, initial]

\(ndsol1\  = \ 
      NDSolve[eq2 //. {b -> 0.05, a -> \(-1\)}

          		, {q[t], p[t]}, {t, 0, 30}] // Flatten\)

\(q1[t_] = q[t] /. ndsol1; 
p1[t_] = p[t] /. ndsol1; \)

\(ndsol2\  = \ 
      NDSolve[eq2 //. {b -> 0.05, a -> \(+4\)}

          		, {q[t], p[t]}, {t, 0, 30}] // Flatten\)

\(q2[t_] = q[t] /. ndsol2; 
p2[t_] = p[t] /. ndsol2; \)

Part (c)  Plot phasee trajectories for the double-well and single-well behavior.

\(\(Show[
      GraphicsArray[{
          Plot[q1[t], {t, 0, 30}
	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
            DisplayFunction -> Identity]

          			\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
          Plot[p1[t], {t, 0, 30}
	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
            DisplayFunction -> Identity]}]

      	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
      DisplayFunction -> $DisplayFunction]; \)\)

\(\(Show[
      GraphicsArray[{
          Plot[q2[t], {t, 0, 30}
	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
            DisplayFunction -> Identity]

          			\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
          Plot[p2[t], {t, 0, 30}
	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
            DisplayFunction -> Identity]}]

      	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
      DisplayFunction -> $DisplayFunction]; \)\)

doublePlot[q1[t], p1[t], {t, 0, 30}];

doublePlot[q2[t], p2[t], {t, 0, 30}];

Part (d)  Animate one of the solutions..

\(frames = 40; 
oscMovie[q2, 30, frames, {\(-15\), 15}]; \)

Problem 2:  Forced Duffing Oscillator for Double-well Potential

Needs["Calculus`FourierTransform`"]

Clear["Global`*`"];

Part (a)  Numerically solve the equations for Q0=0 and Q0=0.1.  Plot q[t] and p[t].

eq1 = 
      q''[t] == 
        a q[t]  -  b q[t]^3  -  γ q'[t]  +  
          Q0 Cos[ω t]

eq2 = eq1 //. {q''[t] -> p'[t], q'[t] -> p[t]}

eq3 = q'[t] == p[t]

eq4 = {eq2, eq3}

initial = {p[0] == 0.001, q[0] == 0}

values = {a -> 0.4, b -> 0.5, ω -> 1/8, γ -> 0.02}

\(ndsol1 = 
	
      NDSolve[\(Join[eq4, initial] //. values\) //. {Q0 -> 0}

          		, {q[t], p[t]}, {t, 0, 120}
	, MaxSteps -> 700] // Flatten
        \)

\(q1[t_] = q[t] /. ndsol1; 
p1[t_] = p[t] /. ndsol1; \)

\(ndsol2 = 
	
      NDSolve[\(Join[eq4, initial] //. values\) //. {Q0 -> 0.1}

          		, {q[t], p[t]}, {t, 0, 120}
	, MaxSteps -> 700] // Flatten
        \)

\(q2[t_] = q[t] /. ndsol2; 
p2[t_] = p[t] /. ndsol2; \)

\(\(Show[
      GraphicsArray[{
          Plot[q1[t], {t, 0, 120}
	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
            DisplayFunction -> Identity]

          			\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
          Plot[p1[t], {t, 0, 120}
	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
            DisplayFunction -> Identity]}]

      	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
      DisplayFunction -> $DisplayFunction]; \)\)

\(\(Show[
      GraphicsArray[{
          Plot[q2[t], {t, 0, 120}
	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
            DisplayFunction -> Identity]

          			\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
          Plot[p2[t], {t, 0, 120}
	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
            DisplayFunction -> Identity]}]

      	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
      DisplayFunction -> $DisplayFunction]; \)\)

Part (b)  Graph the phase trajectory and spectrum for the
       two cases

doublePlot[q1[t], p1[t], {t, 0, 120}];

doublePlot[q2[t], p2[t], {t, 0, 120}];

\(\(plot1 = 
      Show[fftSpectrum[q1[t], t, 120]
	, PlotRange -> {{0, 50}, All}

        	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , DisplayFunction -> Identity]; 
    \)\)

\(\(plot2 = 
      Show[fftSpectrum[q2[t], t, 120]
	, PlotRange -> {{0, 50}, All}

        	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , DisplayFunction -> Identity]; 
    \)\)

\(\(Show[GraphicsArray[{plot1, plot2}]

      	\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ , 
      DisplayFunction -> $DisplayFunction]; \)\)