Method of averaging

• There are two approximations in this what is called first approximation estimate: reduction to the average of the vector field and negligence of ${\displaystyle {\mathcal {O}}(\varepsilon ^{2})}$ terms.
• Uniformity with respect to the initial condition ${\displaystyle x_{0}}$: if we vary ${\displaystyle x_{0}}$ this affects the estimation of ${\displaystyle L}$ and ${\displaystyle c}$. The proof and discussion of this can be found in J. Murdock's book.[3]
• Reduction of regularity: there is a more general form of this theorem which requires only ${\displaystyle f^{1}}$ to be Lipschitz and ${\displaystyle f^{[2]}}$ continuous. It is a more recent proof and can be seen in Sanders et al..[2] The theorem statement presented here is due to the proof framework proposed by Krylov-Bogoliubov which is based on an introduction of a near-identity transformation. The advantage of this method is the extension to more general settings such as infinite-dimensional systems - partial differential equation or delay differential equations.
• J. Hale presents generalizations to almost periodic vector-fields.[4]

Krylov-Bogoliubov realized that the slow dynamics of the system determines the leading order of the asymptotic solution.

In order to proof it, they proposed a near-identity transformation, which turned out to be a change of coordinates with its own time-scale transforming the original system to the averaged one.

1. Determination of a near-identity transformation: the smooth mapping ${\displaystyle y\mapsto U(y,t,\varepsilon )=y+\varepsilon u^{[1]}(y,t,\varepsilon )}$where ${\displaystyle u^{[1]}}$ is assumed to be regular enough and ${\displaystyle T}$ periodic. The proposed change of coordinates is given by ${\displaystyle x=U(y,t,\varepsilon )}$.
2. Choose an appropriate ${\displaystyle u^{[1]}}$solving the homological equation of the averaging theory: ${\displaystyle {\frac {\partial u^{[1]}}{\partial t}}=f^{1}(y,t)-{\bar {f}}^{1}(y)}$.
3. Change of coordinates carries the original system to ${\displaystyle {\dot {y}}=\varepsilon {\bar {f}}^{1}(y)+\varepsilon ^{2}f_{*}^{[2]}(y,t,\varepsilon ).}$
4. Estimation of error due to truncation and comparison to the original variable.

• The same can be done to time-dependent linear parts. Although the fundamental solution may be non-trivial to write down explicitly, the procedure is similar. See Sanders et al. [2] for further details.
• If the eigenvalues of ${\displaystyle A}$ are not all purely imaginary this is called hyperbolicity condition. For this occasion, the perturbation equation may present some serious problems even whether ${\displaystyle g}$ is bounded, since the solution grows exponentially fast.[2] However, qualitatively, we may be able to know the asymptotic solution, such as Hartman-Grobman results and more.[1]
• Occasionally, polar coordinates may yield standard forms that are simpler to analyze. Consider ${\displaystyle z_{1}=r\sin(t-\phi )~{\text{and}}~z_{2}=r\cos(t-\phi )}$, which determines the initial condition ${\displaystyle (r(0),\phi (0))}$ and the system
  $${\displaystyle {\begin{bmatrix}{\dot {r}}\\{\dot {\phi }}\end{bmatrix}}=\varepsilon {\begin{bmatrix}\cos(t-\phi )g(r\sin(t-\phi ),r\cos(t-\phi ),t)\\{\frac {1}{r}}\sin(t-\phi )g(r\sin(t-\phi ),r\cos(t-\phi ),t)\end{bmatrix}}.}$$

  

If ${\displaystyle g\in C^{1}}$we may apply averaging so long as a neighborhood of the origin is excluded (since the polar coordinates fail):

$${\displaystyle {\begin{array}{rcl}{\bar {f}}_{1}^{1}(r)&=&\displaystyle {\frac {1}{2\pi }}\int _{0}^{2\pi }\cos(s-\phi )g(r\sin(s-\phi ),r\cos(s-\phi ),s)ds\\[4pt]{\bar {f}}_{2}^{1}(r)&=&\displaystyle {\frac {1}{2\pi r}}\int _{0}^{2\pi }\sin(s-\phi )g(r\sin(s-\phi ),r\cos(s-\phi ),s)ds,\end{array}}}$$

where the averaged system is

$${\displaystyle {\begin{array}{lcr}{\dot {\bar {r}}}=\varepsilon {\bar {f}}_{1}^{1}({\bar {r}})\\{\dot {\bar {\phi }}}=\varepsilon {\bar {f}}_{2}^{1}({\bar {r}}).\end{array}}}$$

Figure 2: A simple harmonic oscillator with small periodic damping term given by ${\displaystyle {\ddot {z}}+4\varepsilon \cos ^{2}{(t)}{\dot {z}}+z=0,~z(0)=0,~{\dot {z}}(0)=1;~\varepsilon =0.05}$.The numerical simulation of the original equation (blue solid line) is compared with averaging system (orange dashed line) and the crude averaged system (green dash-dotted line). The left plot displays the solution evolved in time and the right plot represents on the phase space. We note that the crude averaging disagrees with the expected solution.

The method contains some assumptions and restrictions. These limitations play important role when we average the original equation which is not into the standard form, and we can discuss counterexample of it. The following example in order to discourage this hurried averaging:[2]

$${\displaystyle {\ddot {z}}+4\varepsilon \cos ^{2}{(t)}{\dot {z}}+z=0,\qquad z(0)=0,\quad {\dot {z}}(0)=1,}$$

where we put ${\displaystyle g(z,{\dot {z}},t)=-4\cos ^{2}(t){\dot {z}}}$following the previous notation.

This systems corresponds to a damped harmonic oscillator where the damping term oscillates between ${\displaystyle 0}$ and ${\displaystyle 4\varepsilon }$. Averaging the friction term over one cycle of ${\displaystyle 2\pi }$ yields the equation:

$${\displaystyle {\ddot {\bar {z}}}+2\varepsilon {\dot {\bar {z}}}+{\bar {z}}=0,\qquad {\bar {z}}(0)=0,\quad {\dot {\bar {z}}}(0)=1.}$$

The solution is

$${\displaystyle {\bar {z}}(t)={\frac {1}{(1-\varepsilon ^{2})^{\frac {1}{2}}}}e^{-\varepsilon t}\sin {((1-\varepsilon ^{2})^{\frac {1}{2}}t)}.}$$

which the convergence rate to the origin is ${\displaystyle \varepsilon }$. The averaged system obtained from the standard form yields:

$${\displaystyle {\begin{array}{lcr}{\dot {\bar {r}}}=-{\frac {1}{2}}\varepsilon {\bar {r}}(2+\cos(2{\bar {\phi }})),~{\bar {r}}(0)=1\\{\dot {\bar {\phi }}}={\frac {1}{2}}\varepsilon \sin(2{\bar {\phi }}),~{\bar {\phi }}(0)=0,\end{array}}}$$

which in the rectangular coordinate shows explicitly that indeed the rate of convergence to the origin is ${\textstyle {\frac {3}{2}}\varepsilon }$ differing from the previous crude averaged system:

$${\displaystyle y(t)=e^{-{\frac {3}{2}}\varepsilon t}\sin {t}}$$

Figure 3: Phase space of a Van der Pol oscillator with ${\displaystyle \varepsilon =0.1}$. The stable limit cycle (orange solid line) in the system is captured correctly by the qualitative analysis of the averaged system. For two different initial conditions ( black dots ) we observe the trajectories.(dashed blue line) converging to the periodic orbit.

Van der Pol was concerned with obtaining approximate solutions for equations of the type

$${\displaystyle {\ddot {z}}+\varepsilon (1-z^{2}){\dot {z}}+z=0,}$$

where ${\displaystyle g(z,{\dot {z}},t)=(1-z^{2}){\dot {z}}}$ following the previous notation. This system is often called the Van der Pol oscillator. Applying periodic averaging to this nonlinear oscillator provides qualitative knowledge of the phase space without solving the system explicitly.

The averaged system is

$${\displaystyle {\begin{array}{lcr}{\dot {\bar {r}}}={\frac {1}{2}}\varepsilon {\bar {r}}(1-{\frac {1}{4}}{\bar {r}}^{2})\\{\dot {\bar {\phi }}}=0,\end{array}}}$$

and we can analyze the fixed points and their stability. There is an unstable fixed point at the origin and a stable limit cycle represented by ${\displaystyle {\bar {r}}=2}$.

The existence of such stable limit-cycle can be stated as a theorem.

Theorem (Existence of a periodic orbit)[5]: If ${\displaystyle p_{0}}$is a hyperbolic fixed point of

$${\displaystyle {\dot {y}}=\varepsilon {\bar {f}}^{1}(y)}$$

Then there exists ${\displaystyle \varepsilon _{0}>0}$such that for all ${\displaystyle 0<\varepsilon <\varepsilon _{0}}$,

$${\displaystyle {\dot {x}}=\varepsilon f^{1}(x,t)+\varepsilon ^{2}f^{[2]}(x,t,\varepsilon )}$$

has a unique hyperbolic periodic orbit ${\displaystyle \gamma _{\varepsilon }(t)=p_{0}+{\mathcal {O}}(\varepsilon )}$ of the same stability type as ${\displaystyle p_{0}}$.

The proof can be found at Guckenheimer and Holmes,[5] Sanders et al. [2] and for the angle case in Chicone.[1]

Figure 4: The plot depicts two fundamental quantities the average technique is based on: the bounded and connected region ${\displaystyle D}$ of the phase space and how long (defined by the constant ${\displaystyle c}$) the averaged solution is valid. For this case, ${\textstyle {\ddot {z}}+z=8\varepsilon \cos {(t)}{\dot {z}}^{2},~z(0)=0,~{\dot {z}}(0)=1;~8\varepsilon ={\frac {2}{15}}}$. Note that both solutions blow up in finite time.  Hence, ${\displaystyle D}$ has been chosen accordingly in order to maintain the boundedness of the solution and the time interval of validity of the approximation is ${\displaystyle 0\leq \varepsilon t<L<{\frac {1}{3}}}$.

The average theorem assumes existence of a connected and bounded region ${\displaystyle D\subset \mathbb {R} ^{n}}$which affects the time interval ${\displaystyle L}$ of the result validity. The following example points it out. Consider the

$${\displaystyle {\ddot {z}}+z=8\varepsilon \cos {(t)}{\dot {z}}^{2},~z(0)=0,~{\dot {z}}(0)=1,}$$

where ${\displaystyle g(z,{\dot {z}},t)=8{\dot {z}}^{2}\cos(t)}$. The averaged system consists of

$${\displaystyle {\begin{array}{lcr}{\dot {\bar {r}}}=3\varepsilon {\bar {r}}^{2}\cos({\bar {\phi }}),~{\bar {r}}(0)=1\\{\dot {\bar {\phi }}}=-\varepsilon {\bar {r}}\sin({\bar {\phi }}),~{\bar {\phi }}(0)=0,\end{array}}}$$

which under this initial condition indicates that the original solution behaves like

$${\displaystyle z(t)={\frac {\sin(t)}{1-3\varepsilon t}}+{\mathcal {O}}(\varepsilon ),}$$

where it holds on a bounded region over ${\displaystyle 0\leq \varepsilon t\leq L<{\frac {1}{3}}}$.

Consider a damped pendulum whose point of suspension is vibrated vertically by a small amplitude, high frequency signal (this is usually known as dithering). The equation of motion for such a pendulum is given by

$${\displaystyle m(l{\ddot {\theta }}-ak\omega ^{2}\sin \omega t\sin \theta )=-mg\sin \theta -k(l{\dot {\theta }}+a\omega \cos \omega t\sin \theta )}$$

where ${\displaystyle a\sin \omega t}$ describes the motion of the suspension point, ${\displaystyle k}$ describes the damping of the pendulum, and ${\displaystyle \theta }$ is the angle made by the pendulum with the vertical.

The phase space form of this equation is given by

$${\displaystyle {\begin{aligned}{\dot {t}}&=1\\{\dot {\theta }}&=p\\{\dot {p}}&={\frac {1}{ml}}(mak\omega ^{2}\sin \omega t\sin \theta -mg\sin \theta -k(lp+a\omega \cos \omega t\sin \theta ))\end{aligned}}}$$

where we have introduced the variable ${\displaystyle p}$ and written the system as an autonomous, first-order system in ${\displaystyle (t,\theta ,p)}$-space.

Suppose that the angular frequency of the vertical vibrations, ${\displaystyle \omega }$, is much greater than the natural frequency of the pendulum, ${\textstyle {\sqrt {g/l}}}$. Suppose also that the amplitude of the vertical vibrations, ${\displaystyle a}$, is much less than the length ${\displaystyle l}$ of the pendulum. The pendulum's trajectory in phase space will trace out a spiral around a curve ${\displaystyle C}$, moving along ${\displaystyle C}$ at the slow rate ${\displaystyle {\sqrt {g/l}}}$ but moving around it at the fast rate ${\displaystyle \omega }$. The radius of the spiral around ${\displaystyle C}$ will be small and proportional to ${\displaystyle a}$. The average behaviour of the trajectory, over a timescale much larger than ${\displaystyle 2\pi /\omega }$, will be to follow the curve ${\displaystyle C}$.