He’s frequency–amplitude formulation with average residuals for nonlinear oscillators

Abstract

In He’s frequency–amplitude formulation, the relationship between frequency and amplitude of a nonlinear oscillator can be obtained through residuals of two trial solutions. Although a high accurate solution can be obtained, this method has some space to be further improved. Here, we show that the calculation of the residuals can be further simplified without loss of accuracy. Duffing oscillator with high nonlinearity is used as an example to show the solution process and accuracy.

Introduction

Nonlinear problems play an important role not only in all areas of physics but also in many other disciplines, such as electronics, chemistry, biomechanics, physiology, textile engineering, and nanotechnology, because fundamental principles of most phenomena in the world are essentially nonlinear. In general, it is very difficult to solve nonlinear problems accurately. After decades of continuous efforts, some excellent methods appeared in open literature, such as the homotopy perturbation methods,^1–11 the variational iteration method,^12–17 the parameter expansion method, ¹⁸ the max–min approach, ¹⁹ the variational approach,^20–23 the Hamiltonian approach, ²⁴ the energy balance method, ²⁵ the stability analysis,^26–29 and the simple formulae for fast prediction of frequency. ³⁰

To solve nonlinear oscillators, Chinese mathematician Prof. Ji-Huan He proposed a simple but effect frequency–amplitude formulation.^31–34 He’s frequency–amplitude formulation was used by many authors with great success.^35–39 In this study, we will suggest an alternative form of the residuals for He’s frequency–amplitude formulation and apply it to high-order Duffing oscillator.

He’s frequency–amplitude formulation

To illustrate the basic idea of He’s frequency–amplitude formulation,^3,31–34 we first consider an algebraic equation in the form

$$f\left(x\right)=0$$

(1)

Letting x₁ and x₂ be two trial solution of equation (1), we obtain two remainders f (x₁) and f (x₂), respectively. According to an ancient Chinese method, its approximate solution can be expressed

$$x=\frac{x_{2}f\left(x_1\right)-x_{1}f\left(x_2\right)}{f\left(x_1\right)-f\left(x_2\right)}$$

(2)

This formula generally gives us an accurate solution, though the calculation procedure is very easy.

Now we consider a generalized nonlinear oscillator with initial conditions in the following form

$$u^{″}+f(u)=0,u(0)=A,u^{\prime }(0)=0$$

(3)

where A is its initial amplitude. According to He’s frequency–amplitude formulation, we use two trial functions $u_1=A\cos \hspace{0.17em}{\omega }_{1}t$ , and $u_2=A\cos \hspace{0.17em}{\omega }_{2}t$ , which are, respectively, the solutions of the following two linear oscillator equations

$$u^{″}+{\mathrm{\omega }}_1^{2}u=0$$

and

$$u^{″}+{\omega }_2^{2}u=0$$

where ω₁ and ω₂ are trial frequencies. The residuals are

$$R_1\text{(}t\text{)}=f\text{(}{u}_1\text{)}-{\omega }_1^{2}A\text{cos}\hspace{0.17em}{\omega }_{1}t$$

(4)

and

$$R_2\text{(}t\text{)}=f\text{(}{u}_2\text{)}-{\omega }_2^{2}A\text{cos}\hspace{0.17em}{\omega }_{2}t$$

(5)

According to classic procedure of He’s frequency–amplitude formulation, the frequency of the nonlinear oscillator equation (3) can be computed approximately by

$${\omega }^2=\frac{{\omega }_1^2{R}_2\left(t_2\right)-{\omega }_2^2{R}_1\left(t_1\right)}{R_2\left(t_2\right)-R_1\left(t_1\right)}$$

(6)

where t₁ and t₂ are location points which determine the phase of the residuals. The residual equations (4) and (5) are sometimes written as

$${\tilde{R}}_1=\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}{R}_1(t)\cos \hspace{0.17em}{\omega }_{1}t\text{d}t}$$

(7)

and

$${\tilde{R}}_2=\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}{R}_2(t)\cos \hspace{0.17em}{\omega }_{2}t\text{d}t}$$

(8)

From these equations, the approximate frequency of equation (3) can also be calculated

$${\omega }_{\mathrm{classic}}^2=\frac{{\omega }_1^2{\tilde{R}}_2-{\omega }_2^2{\tilde{R}}_1}{{\tilde{R}}_2-{\tilde{R}}_1}$$

(9)

Modification of He’s frequency–amplitude formulation

In He’s frequency–amplitude formulation equation (9), the function effects on residuals given in equations (7) and (8) are neglected when the phases ${\omega }_{1}t$ and ${\omega }_{2}t$ are near to $\frac{\pi }{2}$ , while they are overemphasized when the phases are near to 0 in the internal $\text{t}\in {\text{[0}}_{,}\frac{T}{\text{4}}\text{]}$ , thus the weighting factors $\cos {\omega }_{1}t$ and $\cos \hspace{0.17em}{\omega }_{2}t$ given in equations (7) and (8) are not the best ones. In this paper, we will drop the weighting factors $\cos \hspace{0.17em}{\omega }_{1}t$ and $\cos \hspace{0.17em}{\omega }_{2}t$ from the residual equations (7) and (8) to calculate the following residuals

$${\overline{R}}_1=\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}{R}_1(t)\text{d}t}$$

(10)

and

$${\overline{R}}_2=\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}{R}_2(t)\text{d}t}$$

(11)

Equation (9) can be rewritten as

$${\omega }_{\mathrm{Present}}^2=\frac{{\omega }_1^2{\overline{R}}_2-{\omega }_2^2{\overline{R}}_1}{{\overline{R}}_2-{\overline{R}}_1}$$

(12)

Actually, equations (10) and (11) are, respectively, geometric averages of the residuals of equations (4) and (5) in the internal $t\in {\text{[0}}_{,}\frac{T}{\text{4}}\text{]}$ .

Example

To show the effectiveness of our modification, we consider a generalized Duffing oscillator in the form ⁴⁰

$$\ddot{u}+{\alpha }_{0}u+{\alpha }_1{u}^3+{\alpha }_2{u}^5+{\alpha }_3{u}^7+\cdot \cdots \cdot +{\alpha }_n{u}^{2n+1}=0,u(0)=A,\dot{u}(0)=0$$

(13)

Choosing two trial functions $u_1=A\cos \hspace{0.17em}t$ and $u_2=A\cos \hspace{0.17em}2t$ , we have the following two residuals

$$\begin{array}{l}{R}_1(t)=-A\cos \hspace{0.17em}t+{\alpha }_{0}A\cos \hspace{0.17em}t+{\alpha }_1{A}^3\cos 3t+{\alpha }_2{A}^5\cos 5t+\cdots +{\alpha }_n{A}^{2n+1}\cos 2n+1t\\ =({\alpha }_0-1)A\cos \hspace{0.17em}t+{\alpha }_1{A}^3\cos 3t+{\alpha }_2{A}^5\cos 5t+\cdots +{\alpha }_n{A}^{2n+1}\cos 2n+1t\end{array}$$

(14)

and

$$\begin{array}{l}{R}_2\text{(}t\text{)}=-4A\cos \hspace{0.17em}2t+{\alpha }_{0}A\cos \hspace{0.17em}2t+{\alpha }_1{A}^3\cos 32t+{\alpha }_2{A}^5\cos 52t+\cdots +{\alpha }_n{A}^{2n+1}\cos 2n+12t\\ =\text{(}{\alpha }_0-4\text{)}A\cos \hspace{0.17em}2t+{\alpha }_1{A}^3\hspace{0.17em}\cos 32t+{\alpha }_2{A}^5\hspace{0.17em}\cos 52t+\cdots +{\alpha }_n{A}^{2n+1}\hspace{0.17em}\cos 2n+12t\end{array}$$

(15)

where ${\omega }_1=1$ and ${\omega }_2=2$ .

Classic process

Substituting equations (14) and (15) into equations (7) and (8), respectively, the weighted residuals can be calculated

$$\begin{array}{l}{\tilde{R}}_1=\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}{R}_1(t)\cos \hspace{0.17em}t\text{d}t}\\ =\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}\left[({\alpha }_0-1)A{\cos \hspace{0.17em}}^{2}t+{\alpha }_1{A}^3{\cos \hspace{0.17em}}^{4}t+{\alpha }_2{A}^5{\cos \hspace{0.17em}}^{6}t+\cdots +{\alpha }_n{A}^{2n+1}{\cos \hspace{0.17em}}^{2n+2}t\right]\text{d}t}\\ =\frac{1}{2}({\alpha }_0-1)A+\frac{3}{4}\cdot \frac{1}{2}{\alpha }_1{A}^3+\frac{5}{6}\cdot \frac{3}{4}\cdot \frac{1}{2}{\alpha }_2{A}^5+\cdots +\frac{2n+1}{2n+2}\cdot \frac{2n-1}{2n}\cdot \cdots \cdot \frac{5}{6}\cdot \frac{3}{4}\cdot \frac{1}{2}{\alpha }_n{A}^{2n+1}\end{array}$$

and

$$\begin{array}{l}{\tilde{R}}_2=\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}{R}_2\text{(}t\text{)}\cos \hspace{0.17em}2t\text{d}t}\\ =\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}\left[\text{(}{\alpha }_0-4\text{)}A{\cos \hspace{0.17em}}^{2}2t+{\alpha }_1{A}^3{\cos \hspace{0.17em}}^{4}2t+{\alpha }_2{A}^5{\cos \hspace{0.17em}}^{6}2t+\cdots +{\alpha }_n{A}^{2n+1}{\cos \hspace{0.17em}}^{2n+2}2t\right]}\text{d}t\\ =\frac{1}{2}\text{(}{\alpha }_0-4\text{)}A+\frac{3}{4}\cdot \frac{1}{2}{\alpha }_1{A}^3+\frac{5}{6}\cdot \frac{3}{4}\cdot \frac{1}{2}{\alpha }_2{A}^5+\cdots +\frac{2n+1}{2n+2}\cdot \frac{2n-1}{2n}\cdot \cdots \cdot \frac{5}{6}\cdot \frac{3}{4}\cdot \frac{1}{2}{\alpha }_n{A}^{2n+1}\end{array}$$

We get the frequency of the nonlinear oscillator, which is

$$\begin{array}{l}{\omega }_{\mathrm{Classic}}^2=\frac{{\tilde{R}}_2-{\omega }_2^2{\tilde{R}}_1}{{\tilde{R}}_2-{\tilde{R}}_1}\\ ={\alpha }_0+\frac{3}{4}{\alpha }_1{A}^2+\frac{5}{6}\cdot \frac{3}{4}{\alpha }_2{A}^4+\cdots +\frac{2n}{2n+1}\cdot \frac{2n-2}{2n-1}\cdot \cdots \cdot \frac{2}{3}{\alpha }_n{A}^{2n}\end{array}$$

(16)

If $n=2$ , the frequency–amplitude relationship is obtained as

$${\hspace{0.17em}}_{\omega }^{\mathrm{Classic}}=\sqrt{{\alpha }_0+\frac{3}{4}{\alpha }_1{A}^2+\frac{5}{8}{\alpha }_2{A}^4}$$

(17)

if n = 3

$${\hspace{0.17em}}_{\omega }^{\mathrm{Classic}}=\sqrt{{\alpha }_0+\frac{3}{4}{\alpha }_1{A}^2+\frac{5}{8}{\alpha }_2{A}^4+\frac{35}{64}{\alpha }_3{A}^6}$$

(18)

Present process

Now we substitute equations (14) and (15) into equations (10) and (11), respectively, the geometric average of the residuals can be expressed as

$$\begin{array}{l}{\overline{R}}_1=\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}{R}_1\text{d}t}\\ =\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}\left[\left({\alpha }_0-1\right)A\cos \hspace{0.17em}t+{\alpha }_1{A}^3{\cos }^{3}t+{\alpha }_2{A}^5{\cos }^{5}t+\cdots +{\alpha }_n{A}^{2n+1}{\cos }^{2n+1}t\right]}\text{d}t\\ =\frac{2}{\pi }\left[\left({\alpha }_0-1\right)A+\frac{2}{3}{\alpha }_1{A}^3+\frac{4}{5}\cdot \frac{2}{3}{\alpha }_2{A}^5+\cdots +\frac{2n}{2n+1}\cdot \frac{2n-2}{2n-1}\cdot \cdots \cdot \frac{4}{5}\cdot \frac{2}{3}{\alpha }_n{A}^{2n+1}\right]\end{array}$$

and

$$\begin{array}{l}{\overline{R}}_2=\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}{R}_2\text{d}t}\\ =\frac{4}{T}{\displaystyle {\int }_0^{\frac{T}{4}}\left[\left({\alpha }_0-4\right)A\cos \hspace{0.17em}2t+{\alpha }_1{A}^3{\cos }^{3}2t+{\alpha }_2{A}^5{\cos }^{5}2t+\cdots +{\alpha }_n{A}^{2n+1}{\cos }^{2n+1}2t\right]}\text{d}t\\ =\frac{2}{\pi }\left[\left({\alpha }_0-4\right)A+\frac{2}{3}{\alpha }_1{A}^3+\frac{4}{5}\cdot \frac{2}{3}{\alpha }_2{A}^5+\cdots +\frac{2n}{2n+1}\cdot \frac{2n-2}{2n-1}\cdot \cdots \cdot \frac{4}{5}\cdot \frac{2}{3}{\alpha }_n{A}^{2n+1}\right]\end{array}$$

therefore

$$\begin{array}{l}{\omega }_{\mathrm{Present}}^2=\frac{{\overline{R}}_2-4{\overline{R}}_1}{{\overline{R}}_2-{\overline{R}}_1}\\ ={\alpha }_0+\frac{2}{3}{A}^3{\alpha }_1+\frac{4}{5}\cdot \frac{2}{3}{A}^5{\alpha }_2+\cdots +\frac{2n}{2n+1}\cdot \frac{2n-2}{2n-1}\cdot \cdots \cdot \frac{2}{3}{A}^{2n+1}{\alpha }_n\end{array}$$

(19)

If n = 2, the frequency–amplitude relationship obtained as

$${\omega }_{\mathrm{Present}}=\sqrt{{\alpha }_0+\frac{2}{3}{A}^2{\alpha }_1+\frac{8}{15}{A}^4{\alpha }_2}$$

(20)

if n = 3

$${\hspace{0.17em}}_{\omega }^{\mathrm{Presnet}}=\sqrt{{\alpha }_0+\frac{2}{3}{A}^2{\alpha }_1+\frac{8}{15}{A}^4{\alpha }_2+\frac{16}{35}{A}^6{\alpha }_3}$$

(21)

Error analysis

For large values of α₂ or α₃, the relative errors of the approximate frequencies of equations (17), (18), (20), and (21) are, respectively, as follows

$$\underset{{\alpha }_2\to \infty }{\lim }\frac{\left|{\omega }_{\mathrm{Exact}}-{\omega }_{\mathrm{Classic}}\right|}{{\omega }_{\mathrm{Exact}}}=0.0585$$

(22)

$$\underset{{\alpha }_3\to \infty }{\lim }\frac{\left|{\omega }_{\mathrm{Exact}}-{\omega }_{\mathrm{Classic}}\right|}{{\omega }_{\mathrm{Exact}}}=0.0956$$

(23)

$$\underset{{\alpha }_2\to \infty }{\lim }\frac{\left|{\omega }_{\mathrm{Exact}}-{\omega }_{\mathrm{Present}}\right|}{{\omega }_{\mathrm{Exact}}}=0.0221$$

(24)

$$\underset{{\alpha }_3\to \infty }{\lim }\frac{\left|{\omega }_{\mathrm{Exact}}-{\omega }_{\mathrm{Present}}\right|}{{\omega }_{\mathrm{Exact}}}=0.0018$$

(25)

Comparison of the approximate frequencies (ω_Classic and ω_Present) obtained by equations (18) and (21) with the exact frequencies ω_Exact is made in Table 1 for different values of α₁, α₂, and α₃, where the exact frequencies ω_Exact are determined by equation (37) in the Appendix. All present solutions of ω_Present are in good agreement with the exact ones. Other comparisons between the approximate solutions and the exact ones can also be found in Figures 1 and 2, in which the present solution is overlapping with the exact one for A = α₀ =1, α₁ = α₂ = 0, α₃ = 100. From the table and the figures, we could see that the errors of the present solutions ω_Present are decreasing when the values of α₂ and α₃ increase.

Figure 1.

Comparison of present and classic solution with exact one for A = α₀ = α₁ =1, α₂ = α₃ = 0, the present solution gives an exact match with the exact one.

Figure 2.

Comparison of present and classic solution with exact one for A = α₀ =1, α₁ = α₂ = 0, α₃ = 100, the present solution gives an exact match with the exact one.

Table 1.

Relative errors of solutions from present and classic process (A = α₀ = 1).

α 1	α 2	α 3	ω Exact	Present process		Classic process
				ω Present	Relative error	ω Classic	Relative error
1	0	0	1.3178	1.2910	0.0203	1.3229	0.0038
10	0	0	2.8667	2.7689	0.0341	2.9155	0.0170
100	0	0	8.5335	8.2260	0.0360	8.7178	0.0216
1000	0	0	26.804	25.839	0.0360	27.404	0.0224
0	1	0	1.2647	1.2383	0.0209	1.2748	0.0080
0	10	0	2.5837	2.5166	0.0260	2.6926	0.0421
0	100	0	7.5428	7.3711	0.0227	7.9687	0.0565
0	1000	0	23.639	23.116	0.0221	25.020	0.0584
0	0	1	1.2305	1.2071	0.0190	1.2437	0.0107
0	0	10	2.3880	2.3604	0.0116	2.5434	0.0651
0	0	100	6.8362	6.8348	0.0002	7.4624	0.0916
0	0	1000	21.371	21.404	0.0015	23.407	0.0953
1	1	1	1.6753	1.6301	0.0270	1.7094	0.0204
10	10	10	3.4660	3.4017	0.0186	3.6700	0.0589
100	100	100	10.169	10.036	0.0131	10.906	0.0725
1000	1000	1000	31.894	31.498	0.0124	34.258	0.0741

Conclusion

He’s frequency–amplitude formulation modified with geometric averages of residuals can produce quite accurate solutions, especially for high-order Duffing oscillators.