Preface


This section discusses a simplified version of the Adomian decomposition method first concept of which was proposed by Randolf Rach in 1989 that was crystallized later in a paper published with his colleagues G. Adomian and R.E. Meyers. That is way this technique is frequently referred to as the Rach--Adomian--Meyers modified decomposition method (MDM for short). Initially, this method was applied to power series expansions, which was based on the nonlinear transformation of series by the Adomian--Rach Theorem. Similar to the Runge--Kutta methods, the MDM can be implemented in numerical integration of differential equations by one-step methods. In case of polynomials or power series, it shows the advantage in speed and accuracy of calculations when at each step the Adomian decomposition method allows one to perform explicit evaluations.

Return to computing page for the first course APMA0330
Return to computing page for the second course APMA0340
Return to Mathematica tutorial for the second course APMA0340
Return to the main page for the course APMA0330
Return to the main page for the course APMA0340
Return to Part IV of the course APMA0330

Forced Vibrations


 

See http://www.sharetechnote.com/html/DE_Modeling_Example_SpringMass.html#SingleSpring_SimpleHarmonic_Vert_Damp_Force http://www.sharetechnote.com/html/DE_Modeling_Example_SpringMass.html#Single_Spring_with_ExternalForce

This section presents the situation in which a periodic external force is applied to a spring-mass system.

Example: Suppose that the motion of a certain spring-mass system satisfies the differential equation

\[ 4\,\ddot{y} +4\,\dot{y} + 10\,y(t) = 25\,\sin \left( 2\,t \right) \]
and the initial conditions
\[ y(0) =3, \qquad \dot{y} (0) =-1. \]
To find its solution, we first write the characteristic equation \( 4\lambda^2 + 4\lambda +10 =0 . \) Since its nulls are \( \lambda = -\frac{1}{2} \pm {\bf j}\,\frac{3}{2} , \) the general solution of the corresponding homogeneous equation \( 4\,\ddot{y} +4\,\dot{y} + 10\,y(t) = 0 \) is
\[ y_h (t) = c_1 e^{-t/2} \cos \left( \frac{3\,t}{2} \right) + c_2 e^{-t/2} \sin \left( \frac{3\,t}{2} \right) , \]
where c1 and c2 are arbitrary constants. According to the method of undeterming coefficients, we seek a particular solution of the given driven equation in the form
\[ y_p (t) = A\, \cos 2t + B\,\sin 2t , \]
where A and B are found by substituting yp into the given nonhomogeneous equaion. We have \( \dot{y}_p = -2A\,\sin 2t + 2B\,\cos 2t \) and \( \ddot{y}_p = -4A\,\cos 2t - 4B\,\sin 2t = -4\,y_p (t) . \) Thus, upon substituting into the given driven equation, we obtain
\[ \left[ 4\texttt{D}^2 + 4\texttt{D} + 10 \right] y_p = -4\,y_p + 4 \left( -2A\,\sin 2t + 2B\,\cos 2t \right) + 10\, y_p = 25\,\sin 2t , \]
which leads to
\[ \left( 6A +8B \right) \cos 2t + \left( 6B -8A \right) \sin 2t = 25\,\sin 2t . \]
Consequently, A and B must satisfy the equations
\[ 6A +8B =0, \qquad 6B -8A =25 . \]
Using Mathematica, we find
Solve[{6*A + 8*B == 0, 6*B - 8*A == 25}, {A, B}]
Out[1]= {{A -> -2, B -> 3/2}}
Therefore, the particular solution is
\[ y_p (t) = -2\,\cos 2t + \frac{3}{2}\,\sin 2t , \]
and the general solution becomes
\[ y (t) = y_h (t) + y_p (t) = c_1 e^{-t/2} \cos \left( \frac{3\,t}{2} \right) + c_2 e^{-t/2} \sin \left( \frac{3\,t}{2} \right) -2\,\cos 2t + \frac{3}{2}\,\sin 2t . \]
The remaining constants c1 and c2 are determined by the initial conditions that yields
\[ y (0) = c_1 -2 =3, \qquad \dot{y}(0) = - \frac{1}{2}\,c_1 + \frac{3}{2}\, c_2 + 3 = -1. \]
So c1 = 5 and c2 = 3. Then we finally obtain the solution of the given initial value problem:
\[ y(t) = -2\,\cos 2t - \\frac{3}{2}\,\sin 2t + e^{-t/2} \left[ 5\,\cos \left( \frac{3\,t}{2} \right) + 3\, \sin \left( \frac{3\,t}{2} \right) \right] . \]
Finally, we check with Mathematica
DSolve[{4*y''[t] + 4*y'[t] + 10*y[t] == 25*Sin[2*t], y[0] == 3, y'[0] == -1}, y[t], t]
Out[2]= {{y[t] -> -2 Cos[2 t] - 3 Cos[t] Sin[t] + E^(-t/2) (5 Cos[(3 t)/2] + 3 Sin[(3 t)/2])}}

Since our solution required many steps, we show akternative way to find the solution using the Laplace transform. Upon application of the Laplace transformation to the given inital value problem, we reduce it to an algebraic problem:

\[ \left( 4\lambda^2 + 4\lambda + 10 \right) y^L = 8 + 12\lambda + \frac{50}{\lambda^2 + 4} . \]
Therefore, the Laplace transform yL of the (unnown yet) solution is
\[ y^L (\lambda ) = \frac{8 + 12\lambda}{4\lambda^2 + 4\lambda + 10} + \frac{50}{\left(\lambda^2 + 4 \right) \left( 4\lambda^2 + 4\lambda + 10 \right)} . \]
Application of the inverse Laplace transform yields the solution:
InverseLaplaceTransform[(8 + 12*s)/(4*s^2 + 4*s + 10), s, t]
FullSimplify[ComplexExpand[%]]
Out[4]= 1/3 E^(-t/2) (9 Cos[(3 t)/2] + Sin[(3 t)/2])
and
InverseLaplaceTransform[(50)/(4*s^2 + 4*s + 10)/(s^2 + 4), s, t]
FullSimplify[ComplexExpand[%]]
Out[6]= -2 Cos[2 t] - 3 Cos[t] Sin[t] + 2/3 E^(-t/2) (3 Cos[(3 t)/2] + 4 Sin[(3 t)/2])
Adding these two expression, we get
2 Cos[2 t] - 3 Cos[t] Sin[t] + 2/3 E^(-t/2) (3 Cos[(3 t)/2] + 4 Sin[(3 t)/2]) + 1/3 E^(-t/2) (9 Cos[(3 t)/2] + Sin[(3 t)/2])
FullSimplify[%]
Out[7]= -2 Cos[2 t] - 3 Cos[t] Sin[t] + E^(-t/2) (5 Cos[(3 t)/2] + 3 Sin[(3 t)/2])
which is exactly the same expression as we already obtained.

It is important to note that the solution consists of two distinct parts. the terms containing the exponential factor \( e^{-t/2} \) rapidly approach zero. It is customary to call these terms transient. the remaining terms involve only sines and cosines, so they represent an oscillation that last indefinitly. We refer to them as a steady state. The transient part comes from the solution of the corresponding homogeneous equation \( 4\,\ddot{y} +4\,\dot{y} + 10\,y(t) = 0 \) subject to the given inital conditions. It is presented on the graph by red curve. The steady state is the particular solution of the driven equation \( 4\,\ddot{y} +4\,\dot{y} + 10\,y(t) = 25\,\sin \left( 2\,t \right) . \) Shortly, the transient solution vanishes and the full solution is essentially indistinguishable from the steady state.

y[t_] = -2 Cos[2 t] - 3 Cos[t] Sin[t] + E^(-t/2) (5 Cos[(3 t)/2] + 3 Sin[(3 t)/2])
yp[t_] = E^(-t/2) (5 Cos[(3 t)/2] + 3 Sin[(3 t)/2])
yh[t_] = -2 Cos[2 t] - 3 Cos[t] Sin[t]
Plot[{y[t], yp[t], yh[t]}, {t, 0, 15}, PlotStyle -> {Thick, Red, Dashed}]

 

 

  1. Cveticanin, L.,Strongly Nonlinear Oscillators–Analytical Solutions, Undergraduate Lecture Notes in Physics, Springer, 2014.
  2. Cveticanin, L., Luburic, U.K., Mester, G., Periodic Motion in an Excited and Damped Cubic Nonlinear Oscillator, Mathematical Problems in Engineering, 2018, Volume 2018, Article ID 3841926, 6 pages; https://doi.org/10.1155/2018/3841926
  3. Cveticanin, L., Mester, G., Biro, I., Parameter influence on the harmonically excited Duffing oscillator, Acta Polytechnica Hungarica, vol. 11, no. 5, pp. 145–160, 2014.
  4. Holmes, P., A nonlinear oscillator with a strange attractor, Philosophical Transactions of the Royal Society A: mathematical, Physical, and Engineering Sciences, 1997, Volume 292, Issue 1394, pp. https://doi.org/10.1098/rsta.1979.0068
  5. Hsu, C.S., On the application of elliptic functions in non-linear forced oscillations, Quarterly of Applied Mathematics, vol. 17, pp. 393–407, 1959/1960.

 

Return to Mathematica page
Return to the main page (APMA0330)
Return to the Part 1 (Plotting)
Return to the Part 2 (First Order ODEs)
Return to the Part 3 (Numerical Methods)
Return to the Part 4 (Second and Higher Order ODEs)
Return to the Part 5 (Series and Recurrences)
Return to the Part 6 (Laplace Transform)
Return to the Part 7 (Boundary Value Problems)