Paradoxical Stabilization of Forced Oscillations by Strong Nonlinear Friction

TL;DR

This paper demonstrates that strong nonlinear friction can paradoxically stabilize forced oscillations by counteracting the ponderomotive force, leading to a stable oscillation state.

Contribution

It reveals a novel stabilization mechanism where nonlinear friction stabilizes oscillations against inhomogeneity-induced forces.

Findings

Strong nonlinear friction creates a stabilizing tug opposite to the ponderomotive force.

Friction can stabilize oscillations near amplitude maxima.

The stabilization occurs when friction-induced force exceeds the ponderomotive force.

Abstract

In a dissipative dynamic system driven by an oscillating force, a strong nonlinear highly oscillatory friction force can create a quasi-steady tug, which is always directed opposite to the ponderomotive force induced due to a spatial inhomogeneity of oscillations. When the friction-induced tug exceeds the ponderomotive force, the friction stabilizes the system oscillations near the maxima of the oscillation spatial amplitude.

Equations

$$\ddot{x}+\mathcal{K}(\mathcal{F})\dot{x}=\mathcal{F},\mathcal{K}=\nu {\mathcal{F}}^{2n}.$$

$$\mathcal{F}(x,t)=f_1(x)\cos \omega t+f_2(x)\sin \omega t.$$

$$x(t)=X(t)+\xi (t)$$

$$⟨\xi ⟩=(\omega /2\pi ){\int }_0^{2\pi /\omega }\xi (t)dt=0.$$

$$\mathcal{F}\approx F+\xi {\partial }_{X}F,F=\mathcal{F}(X,t);$$

$$\displaystyle\mathcal{K}\approx K+\xi\partial_{X}K,\quad K=\mathcal{K}\big{(}\mathcal{F}(X,t)\big{)},$$

$$\ddot{X}+\ddot{\xi }+K\dot{X}+K\dot{\xi }+\xi \dot{X}{\partial }_{X}K+\xi \dot{\xi }{\partial }_{X}K=F+\xi {\partial }_{X}F,$$

$$\kappa =⟨K⟩>0.$$

$$\ddot{\xi }+K\dot{\xi }=F.$$

$$\xi =e^{-\kappa t-\Delta (t)}{\int }_0^t{e}^{\kappa \tau +\Delta (\tau )}F(\tau )d\tau ,$$

$$\Delta (t)={\int }_0^t[K(\tau )-\kappa ]d\tau .$$

$$\xi =\frac{(\kappa f_1-\omega f_2)\sin \omega t-(\kappa f_1+\omega f_2)\cos \omega t}{\omega ({\kappa }^2+{\omega }^2)}$$

$$\ddot{X}+(\kappa +⟨\xi {\partial }_{X}K⟩)\dot{X}=⟨\xi {\partial }_{X}F⟩-⟨\xi \dot{\xi }{\partial }_{X}K⟩.$$

$$\ddot{X}$$

$$\kappa =2^{-2n}\left(\genfrac{}{}{0ex}{}{2n}{n}\right)\nu (f_1^2+f_2^2{)}^n,$$

$$\nu >\frac{2^{2n}(n!{)}^2(n+1{)}^{1/2}}{(2n)!(4n^2-n-1{)}^{1/2}}\frac{\omega }{f_1^{2n}}.$$

$$t=\tau /\omega ,x(t)=(f_0/{\omega }^2)y(\tau );$$

$$\sigma =\nu f_0^{2n}/\omega .$$

$$y^{\prime \prime }(\tau )+\sigma {\cos }^{2n}(\tau )y^{\prime }(\tau )=\cos \tau .$$

$$y^{\prime }(\tau )=e^{-\sigma S_n(\tau )}[y^{\prime }(0)-Y_n(0)]+Y_n(\tau ),$$

$$Y_n(\tau )=\frac{e^{-\sigma S_n(\tau )}}{e^{\pi \alpha \sigma }+1}\underset{-\pi /2}{\overset{\pi /2}{\int }}{e}^{\sigma S_n(\zeta +\tau +\frac{\pi }{2})}\sin (\zeta +\tau )d\zeta ,$$

$$S_n(\tau )=\alpha \tau +\sum _{m=1}^n\left(\genfrac{}{}{0ex}{}{2n}{n+m}\right)\frac{\sin (2m\tau )}{2^{2n}m},\alpha =2^{-2n}\left(\genfrac{}{}{0ex}{}{2n}{n}\right).$$

$$A_n=\frac{1}{\cosh (\frac{\pi \alpha \sigma }{2})}\underset{-\pi /2}{\overset{\pi /2}{\int }}{e}^{\sigma S_n(\tau )}\cos (\tau )d\tau .$$

$$Y_1(\tau )=\sum _{m=1}^{\infty }\left[C_m{e}^{i(2m-1)\tau }+C_m^{\ast }{e}^{-i(2m-1)\tau }\right],$$

$$C_1=\sum _{k=-\infty }^{\infty }\frac{(-1{)}^{k+1}}{4k+2-i\sigma }\left[i{I}_k^2\left(\frac{\sigma }{4}\right)+I_k\left(\frac{\sigma }{4}\right)I_{k+1}\left(\frac{\sigma }{4}\right)\right],$$

$$C_m=i^m\frac{I_{m-\frac{1}{2}-\frac{i\sigma }{4}}\left(\frac{\sigma }{4}\right)}{I_{\frac{3}{2}-\frac{i\sigma }{4}}\left(\frac{\sigma }{4}\right)}\left[\left(2+\frac{4i}{\sigma }\right)C_1+C_1^{\ast }-\frac{2}{\sigma }\right].$$

$$C_1(\sigma )\approx \frac{\sigma /2}{4+{\sigma }^{2-\delta }}-i\left[\frac{\sigma /3-2}{4+{\sigma }^{2-\delta }}+\frac{16+3{\sigma }^2}{(4+{\sigma }^{2-\delta }{)}^2}\right]$$

$$F=f_1\cos (\omega t)+f_2\sin (\omega t),$$

$$K=\nu {\left[f_1\cos (\omega t)+f_2\sin (\omega t)\right]}^{2n},$$

$$\xi =\frac{(\kappa f_1-\omega f_2)\sin \omega t-(\kappa f_1+\omega f_2)\cos \omega t}{\omega ({\kappa }^2+{\omega }^2)},$$

$$\ddot{X}+(\kappa +⟨\xi {\partial }_{X}K⟩)\dot{X}=⟨\xi {\partial }_{X}F⟩-⟨\xi \dot{\xi }{\partial }_{X}K⟩.$$

$$\displaystyle\begin{split}&\kappa=\langle K\rangle=\frac{\omega}{2\pi}\int_{0}^{2\pi/\omega}\nu\left[f_{1}\cos(\omega t)+f_{2}\sin(\omega t)\right]^{2n}dt=\frac{1}{2\pi}\int_{0}^{2\pi}\nu\left[\frac{f_{1}-if_{2}}{2}e^{i\tau}+\frac{f_{1}+if_{2}}{2}e^{-i\tau}\right]^{2n}d\tau=\\ &=\frac{\nu}{2\pi}\sum\limits_{m=0}^{2n}\binom{2n}{m}\left(\frac{f_{1}+if_{2}}{2}\right)^{m}\left(\frac{f_{1}-if_{2}}{2}\right)^{2n-m}{\color[rgb]{1,0,0}\int_{0}^{2\pi}e^{2(n-m)i\tau}d\tau}=\frac{\nu}{2\pi}\binom{2n}{n}\left(\frac{f_{1}+if_{2}}{2}\right)^{n}\left(\frac{f_{1}-if_{2}}{2}\right)^{n}{\color[rgb]{1,0,0}2\pi}.\end{split}$$

$$\kappa =2^{-2n}\left(\genfrac{}{}{0ex}{}{2n}{n}\right)\nu (f_1^2+f_2^2{)}^n.$$

$$⟨\xi {\partial }_{X}K⟩=\frac{\omega }{2\pi }{\int }_0^{2\pi /\omega }\xi (t)\times 2n\nu {\left[f_1\cos (\omega t)+f_2\sin ({\omega }_t)\right]}^{2n-1}\left[{\partial }_X{f}_1\cos (\omega t)+{\partial }_X{f}_2\sin (\omega t)\right]dt=0.$$

$$\begin{array}{r}⟨\xi {\partial }_{X}F⟩=\frac{\omega }{2\pi }{\int }_0^{2\pi /\omega }\frac{(\kappa f_1-\omega f_2)\sin \omega t-(\kappa f_1+\omega f_2)\cos \omega t}{\omega ({\kappa }^2+{\omega }^2)}\left[f_1\cos (\omega t)+f_2\sin ({\omega }_t)\right]dt=\\ =-\frac{f_1{\partial }_X{f}_1+f_2{\partial }_X{f}_2}{2({\kappa }^2+{\omega }^2)}-\frac{\kappa (f_2{\partial }_X{f}_1-f_1{\partial }_X{f}_2)}{2\omega ({\kappa }^2+{\omega }^2)}.\end{array}$$

Full text

I Introduction

In classical mechanics, Kapitza pendulum or Stephenson-Kapitza pendulum Kapitza1 ; Kapitza2 is a statically unstable inverted pendulum whose statically unstable equilibrium position is stabilized by small fast vertical oscillations of the pivot point. This induced stability initially described by A. Stephenson Stephenson has been first explained by P. L. Kapitza Kapitza1 . To find theoretical reason of the induced stability Kapitza separated the pendulum motion into fast and slow oscillations and, by averaging out fast ones, found the effective potential which has minimum at the pendulum upper position, in contrast to a simple pendulum. This approach created a new concept of dynamic stabilization in mechanics Blekhman . Chelomei’s pendulum provides another well known example of dynamically stabilized mechanical system Chelomei . The induced stability can be vindicated beyond the framework of the method of fast and slow motions separation Butikov .

In general, a spatial inhomogeneity of an oscillating driving force creates a quasi-steady ponderomotive force LLM ; Licht-Lieb , directed against the spatial gradient of the driving force. In the Kapitza pendulum, this ponderomotive force acts against gravitation and makes the upper position stable. A dissipation dampens oscillations around the upper equilibrium position thus further stabilizing it. When no other forces present besides the oscillating driving force, as in the case of charged particle dynamics in a standing electromagnetic wave, the ponderomotive force always repels particles from the maxima of the wave spatial amplitude Blekhman .

Here we present a general model of a dissipative dynamic system driven by an oscillating force, where a strong nonlinear highly oscillatory friction creates a quasi-steady tug, which, quite counter-intuitively, is always directed opposite to the ponderomotive force and exceeds it for a sufficiently strong driving force. This leads to a seemingly paradoxical stabilization of the system oscillations near the maxima of the spatial amplitude of the driving force. It differs from the Kapitza pendulum effect in that here the stabilization factor is a nonlinear growth of the friction with the driver force, which creates a tug against the ponderomotive potential.

II Model

We consider a simple one-dimensional model of a forced oscillation with a strong nonlinear friction, given by the equation

[TABLE]

Here the dot denotes differentiation with respect to time; $n$ is a natural number. The friction coefficient $\mathcal{K}$ is a non-negative function of the oscillating driving force $\mathcal{F}$,

[TABLE]

The model is motivated by the dynamics of a charged particle in a strong electromagnetic field, where accelerating particles lose energy and undergo a recoil due to their emission of electromagnetic radiation Rad-Dom . This causes a friction which nonlinearly grows with the electromagnetic field strength.

Following the classical approach of Ref. LLM , we assume that a solution to Eq. (1) can be represented as

[TABLE]

with a slowly varying function $X(t)$ and a fast oscillating small addition $\xi(t)$, $|\xi|\ll|X|$, which has a zero time average,

[TABLE]

Correspondingly, $\langle x\rangle=X$, $\langle\dot{\xi}\rangle=\langle\ddot{\xi}\rangle=0$, $\left<X\right>\approx X(t)$. Substituting (3) into Eq. (1) and expanding the functions $\mathcal{F}$ and $\mathcal{K}$ in powers of $\xi$ as

[TABLE]

we obtain

[TABLE]

where $\partial_{X}$ is the partial derivative with respect to $X$ (the first argument of $\mathcal{F}$). The time derivatives $\ddot{\xi}$ and $\dot{\xi}$ are not small, being proportional to $\omega^{2}$ and $\omega$, respectively. They are assumed to be much greater than $\ddot{X}$ and $\dot{X}$. The friction coefficient $K$ defined in Eq. (6) is not necessarily small; it has a time-independent component

[TABLE]

In Eq. (7), slowly varying and fast oscillating terms should cancel out separately. Neglecting the time derivatives of $X$, in the zeroth order approximation with respect to $\xi$ we find for the fast oscillating term

[TABLE]

Here $X$ as an argument of functions $F$ and $K$ is assumed to be constant. The forced oscillation solution of Eq. (9) can be cast in the form

[TABLE]

The first term in the expansion with respect to $\Delta$ of the dependence given by Eq. (10)

[TABLE]

approximates the first harmonic of the solution.

Averaging Eq. (7) over time and taking into account that $\left<F(X,t)\right>\approx 0$ for nearly constant $X(t)$, we obtain

[TABLE]

Using here the expression for $\xi$, given by Eq. (12), and the definitions for $F$ and $K$ formulated above, we obtain the equation for the slowly varying function $X(t)$ and the average friction coefficient:

[TABLE]

where $\binom{2n}{n}$ is the binomial coefficient. The first term on the r.h.s. of Eq. (14) corresponds to the classical ponderomotive force LLM , modified due to the friction, the last terms represent the tug induced by the friction. In the case of $f_{2}=0$, the third term in the r.h.s. of Eq. (14) vanishes. Then the tug becomes always directed opposite to the ponderomotive force. It can even exceed the latter in magnitude, when

[TABLE]

In this case, trajectories of the system described by the model (1) drift to the local maximum of the driving force spatial amplitude. If the system were non-dissipative, the oscillations near that maximum would be destabilized by the ponderomotive force, so that the corresponding trajectories would drift against the spatial gradient of the driving force and would eventually escape to regions of a lower spatial amplitude of the driving force Blekhman . Sufficiently strong friction makes the trajectories to drift along the spatial gradient of the driving force, provided that these trajectories have already got to a region with sufficiently high driving force. This causes a seemingly paradoxical stabilization of the oscillations near the local maximum of the driving force spatial amplitude.

III Numerical simulations

The effect of the stabilization of oscillations due to strong friction is further demonstrated by the numerical integration of the model equation (1) with the friction coefficient taken in the form $\mathcal{K}=\nu\mathcal{F}(x,t)^{4}$ in two cases. In the first case, shown in Fig. 1 (a,b,c), the driving force amplitude is bell-shaped, $\mathcal{F}(x,t)=f_{0}\exp(-(x/l_{0})^{2})\cos\omega t$, with the width, oscillation frequency, amplitude, and friction factor equal to $l_{0}=10$, $\omega=1$, $f_{0}=3$, and $\nu=0.2$, respectively. In the second case, Fig. 1 (d,e,f), the driving force is spatially periodic, $\mathcal{F}(x,t)=f_{0}\cos^{2}(2\pi x/l_{0})\cos\omega t$, $l_{0}=10$, $\omega=1$, $f_{0}=8$, and $\nu=0.25$.

In the case of bell-shaped driving force shown in Fig. 1 (a,b,c), the trajectories, starting at $t=0$ from locations where the driving force spatial amplitude is relatively weak, exhibit several oscillations and then escape being pushed away by the ponderomotive force. In the region of a high spatial amplitude of the driving force, the friction-induced tug overcomes the ponderomotive force, therefore the trajectories started from this region drift towards the maximum of the driving force spatial amplitude. As they get closer to that maximum, their drift becomes slower and their oscillation amplitude decreases.

In the case of spatially periodic driving force (see Fig. 1 (d,e,f)), the trajectories, initially oscillating near the minima of the driving force spatial amplitude and having enough large oscillation amplitude, reach the region of higher driving force, where the trajectories are caught by the friction-induced tug. Eventually all such trajectories are reduced to oscillations near the maxima of the driving force spatial amplitude.

In both cases the drift is slowing down near the maxima of the driving force spatial amplitude in agreement with Eq. (14), because the gradient $\partial_{X}f$ vanishes at the maximum of $f$.

IV Limit cycle

On a trajectory asymptotically turning into periodic oscillations seen in Fig. 1 (b,c,e,f), the driving force is almost constant, which simplifies theoretical consideration allowing more detail description of the driven oscillations. This case is described by the approximation of $\mathcal{F}(x,t)=f_{0}\cos\omega t$. We change variables to $(\tau,y)$ and introduce the friction parameter $\sigma$:

[TABLE]

Then Eq. (1) becomes

[TABLE]

The general solution can be represented as

[TABLE]

As one can see, any solution at $\tau\rightarrow\infty$ tends to the limit cycle described by the function ${Y}_{n}$.

The amplitude of the limit cycle in terms of the derivative $y^{\prime}$, ${A}_{n}=\max|{Y}_{n}(\tau+\pi)-{Y}_{n}(\tau)|$, is

[TABLE]

For large $\sigma$, it decreases as a negative power of $\sigma$, e.g., ${A}_{1}\approx 1.88\sigma^{-2/3}$. Correspondingly, the oscillation amplitude decreases when the trajectory drifts towards the driving force maximum, as seen in Fig. 1.

In the case of $n=1$, the function describing the limit cycle, Eq. (21), can be represented as a Fourier series in terms of odd harmonics of the driving force

[TABLE]

where the asterix denotes complex conjugation. The sequence $C_{1}$, $C_{m}$, $m\geq 2$, representing the frequency spectrum of ${Y}_{1}$, is expressed in terms of modified Bessel functions of the first kind, $I_{k}$:

[TABLE]

The amplitude of the first harmonic, $C_{1}$, as a function of $\sigma$ can be approximated by

[TABLE]

with $\delta\approx 0.098$. The spectral density of the function ${Y}_{1}(\tau)$, Eq. (24), is $|2C_{m}|^{2}$; it is shown in Fig. 2. As one can see, the frequency spectrum contains high order harmonics according to Eqs. (24-26).

V Conclusion

In conclusion, in contrast to known dynamical destabilization under the action of dissipation known as dissipation-induced instabilities (see review article Destabilization and references therein) we show that a strong nonlinear friction can cause a seemingly paradoxical stabilization of forced oscillations near the maxima of the driving force spatial amplitude. In particular, such a friction occurs in the dynamics of charged particles in ultra-strong electromagnetic fields due to radiation reaction Rad-Dom . The threshold for the described stabilization of forced oscillations corresponds to the criterion of the importance of the radiation reaction effects Survey .

In a standing electromagnetic waves, which can be formed in multiple high power laser configurations Multi-las , the stabilization due to a strong radiation reaction is manifested in an anomalous electron bunching near the electric field maxima Anom-bunch . As shown in Refs. Anom-bunch ; Fedotov ; Attr ; Jirka ; Kirk , electrons can be captured for many laser periods due to radiation friction impeding the ponderomotive force. When radiation reaction dominates, the electron motion in a standing wave evolves to limit cycles and strange attractors Attr ; Jirka . In the case of a circularly polarized standing wave, analytical expressions exist for the limit cycles near the electric field maxima Kirk . A collision of multiple ultra-intense electromagnetic waves creates structurally determinate patterns in the electron phase space Survey ; Vranic ; Gong due to a counterplay of the ponderomotive force and the friction-induced tug.

Although in the present work we have been motivated by the intention to build up the theory of the radiative electron dynamics in the field of extremely high intensity lasers, we believe that the formulated above concept of dissipative stabilization of nonlinear dynamic systems will be useful for applications well beyond the framework of the laser-matter interaction physics MTB ; Marklund , remembering a saying of William Thomson (Lord Kelvin) “I never satisfy myself until I can make a mechanical model of a thing” Kelvin .

Appendix

Here we present mathematical derivations for some formulae shown above. Notations for variables are the same as in the main text. The equation numbering is preserved for those equations which appear in the main text; auxiliary formulae are numbered within sections.

V.1 Equations (14) and (15)

In this section we derive Eqs. (14) and (15) from Eq. (13). We use the following formulae explicitly written or assumed in the main text:

[TABLE]

The time-averaged friction coefficient, defined in Eq. (8) in the main text, is

[TABLE]

Here $\binom{2n}{m}=(2n)!/[m!(2n-m)!]$ is the binomial coefficient. The last integral (marked with red) is zero except the case $m=n$, for which it equals $2\pi$ (therefore the sum contains the only nonzero term, for the index $m=n$). In this way we obtain Eq. (15) from the main text:

[TABLE]

In Eq. (13), the term $\langle\xi\partial_{X}K\rangle$ is zero because it contains only odd harmonics,

[TABLE]

The next averaged term is obtained by a simple integration

[TABLE]

The last averaged terms is

[TABLE]

Here $\delta_{m,n}$ is the Kronecker delta; it equals one for $m=n$ and zero otherwise. In the last line $\nu$ is expressed in terms of $\kappa$ using Eq. (15). Combining the results into Eq. (13) we obtain Eq. (14) from the main text:

[TABLE]

V.2 Equations (20-22)

In this section we derive Eqs. (20), (21) and (22), representing the solution of Eq. (19):

[TABLE]

We seek the solution for the first derivative $y^{\prime}$ in the form

[TABLE]

Substituting this ansatz into Eq. (19), we obtain

[TABLE]

Integrating this equation for $M$ and substituting the result into Eq. (V.2.1), we obtain a general solution of Eq. (19):

[TABLE]

The term in the exponent can be expanded as

[TABLE]

Here we rearranged the sums changing to new indices (marked with red) and used the identity $\binom{2n}{n-m}=\binom{2n}{n+m}$ for ${n\geq m\geq 0}$. In this way we obtain Eq. (22) in the main text:

[TABLE]

Below we use the following properties of the function $S_{n}$:

[TABLE]

The limit cycle is a periodic solution $Y_{n}$ of Eq. (19), therefore for any $\tau$

[TABLE]

Using (V.2.3) we rewrite this expression as

[TABLE]

which gives the formula for $Y_{n}(\tau)$

[TABLE]

In this way we obtain Eq. (21) in the main text:

[TABLE]

Accepting $Y_{n}$ as a particular solution of Eq. (19) and rewriting Eq. (V.2.3), we obtain another form of the general solution of Eq. (19), namely Eq. (20):

[TABLE]

V.3 Equation (23)

In this section we derive Eq. (23). The maximum difference $|{Y}_{n}(\tau+\pi)-{Y}_{n}(\tau)|$ is reached for $\tau=-\pi/2$, therefore

[TABLE]

Thus we obtain Eq. (23) in the main text:

[TABLE]

V.4 Equations (25) and (26)

In this section we derive Eqs. (25) and (26) for the case of $n=1$. The function $Y_{1}$ is a periodic solution of the equation

[TABLE]

It is represented as a Fourier series, which obviously should contain only odd hamonics of the driver:

[TABLE]

where symbol “*” denotes complex conjugation. Substituting Eq. (24) into Eq. (V.4.1) we obtain

[TABLE]

Rearranging the sums and collecting the Fourier coefficients of the terms corresponding to the same harmonics, we easily find the following recursive relations:

[TABLE]

The last expression is the same as a recurrence identity for the modified Bessel function of the first kind math ,

[TABLE]

for $z=\sigma/4$ and $\mu=m-1/2-i\sigma/4$, and $I_{\mu}(z)\propto i^{-m}C_{m}$. The general solution of Eq. (V.4.4) is a linear combination of the modified Bessel functions of the first and second kinds, where one should cancel out an unbounded term:

[TABLE]

The coefficient of proportionality, $U$, is determined using the relation Eq. (V.4.3):

[TABLE]

As a result we obtain Eq. (26) of the main text:

[TABLE]

The formula Eq. (25) for the coefficient $C_{1}$ can be obtained using series representation for the function $Y_{1}$. For $n=1$, we have $\alpha=1/2$, $S_{1}(\tau)=\tau/2+\sin(2\tau)/4$, and Eq. (21) becomes

[TABLE]

Using the generating function for the modified Bessel function of the first kind math ,

[TABLE]

we change exponent terms involving $\sin$ function into series in the following way

[TABLE]

Then the right-hand side of Eq. (V.4.8) is transformed into a double series

[TABLE]

Integrating term-by-term and rearranging the sums by collecting terms involving the same harmonics, it is not difficult to obtain the following formulae

[TABLE]

Using Eq. (V.4.13) and the symmetry property $I_{-k}(z)=I_{k}(z)$ math , it is easy to show that $C_{1-m}=C^{*}_{m}$. For $m=1$ we obtain Eq. (25) from the main text:

[TABLE]

The representation Eq. (V.4.13) is equivalent to Eq. (26) provided that $C_{1}$ is defined by Eq. (25).

Fig. 3 shows a comparison of the numerical solution of Eq. (V.4.1) for $\sigma=20$ with the analytical solution defined by Eqs. (24), (25), and (26), where the sum in Eq. (24) is cut at $m=6$ (which corresponds to the $2m-1=11^{\rm th}$ harmonic) and the sum in Eq. (25) is cut at $k=\pm 10$ (i.e. only terms with the index $k$ from -10 to 10 are taken into account).