Semi-Analytical Closed-Form Solutions of the Ball–Plate Problem

Abstract

Mathematical models and numerical simulations are necessary to understand the dynamical behaviors of complex systems. The aim of this work is to investigate closed-form solutions for the ball–plate problem considering a system derived from an optimal control problem for ball–plate dynamics. The nonlinear properties of ball and plate control system are presented in this work. To semi-analytically solve this system, we explored a second-order nonlinear differential equation. Consequently, we obtained the approximate closed-form solutions by the Optimal Parametric Iteration Method (OPIM) using only one iteration. A comparison between the analytical and corresponding numerical procedures reflects the advantages of the first one. The accordance between the obtained results and the numerical ones highlights that the procedure used is accurate, effective, and good to implement in applications such as sliding mode control to the ball-and-plate problem.

1. Introduction

Chaotic dynamic systems have an application in various sectors, including medicine, mechatronics and secure communication. Therefore, the study of chaotic dynamic systems has real and strong motivation. Jurdjevic [1] has introduced the ball–plate problem as an optimal control problem on the Lie group ${\mathbb{R}}^2\times SO\left(3\right)$. Challenges related to the control and fast dynamic response requiring fast and short sensing and the immediate correction of the selected controller represent the recent development of ball–plate system applications. A mechatronic learning platform of a ball–plate system constitutes a successful tool for training in robotics, automation control methods and applications due to importance of controlling fast and unstable systems [2].

Refs. [3,4,5,6,7] introduce many iteration methods for solving different nonlinear problems. El-Amrani et al. [8] presented an iterative scheme for the numerical analysis modeling of flow and heat transfer in porous media, taking into account the nonlinear viscosity and diffusion coefficient depending on the temperature. Li [9] studied the nonlinear Schrodinger equation with a highly oscillatory potential by means of Picard’s iterative schemes. Hernandez-Veron et al. [10] developed an accurate iterative scheme for solving Wiener–Hopf problems. Hernandez-Veron et al. [11] used Newton-type iterative schemes to analyze the existence of solutions and approximate solutions of Chandrasekhar H-equations. Sharma et al. [12] proposed a new three-step iteration scheme with applications to different problems. Liu et al. [13] applied three novel fifth-order iterative schemes in practical applications. Cordero et al. [14] proposed an iterative scheme to obtain multiple solutions simultaneously. Kapoor et al. [15] obtained series solutions to fractional telegraph equations using an iterative scheme based on Yang’s transform. Ghosh et al. [16] analytically solved some arbitrary-order nonlinear Volterra integro-differential equations involving delay using an iterative scheme. Botchev et al. [17] solved some nonlinear heat-conduction problems by means of an adaptive iterative explicit scheme. Mastryukov [18] presented a different scheme for wave equations using the Laguerre transform. Praks et al. [19] solved the implicit Colebrook equation for flow friction using an iterative scheme, namely the optimal multi-point iterative method.

Other analytical methods, such as the Optimal Parametric Iteration Method (OPIM) [20], the Optimal Homotopy Asymptotic Method (OHAM) [21,22,23], the Optimal Homotopy Perturbation Method (OHPM) [24,25,26] and the modified Optimal Parametric Iteration Method [27], approach nonlinear differential equations. Based on these methods, we performed an analytical study of the ball–plate problem using the OPIM technique. This analytical technique does not depend on the existence of small parameters in equations or initial/boundary conditions.

The dynamic properties of chaotic systems for which the exact solution can no longer be controlled are studied using specific mathematical methods, such as bifurcation diagrams, bifurcation routes, Hopf bifurcation, the existence of heteroclinic orbits or homoclinic orbits, topological horseshoes, Poincare maps, Lyapunov exponent spectra, phase portraits, frequency spectra, equilibria, a dissipative system, and amplitude modulations.

The introduction, challenges, and directions are presented in Section 1, while the ball–plate problem and closed-form solutions are given and discussed in Section 2. Section 3 describes the Optimal Parametric Iteration Method (OPIM) and semi-analytical solutions. Section 4 contains the numerical results. Finally, Section 5 is devoted to some conclusions and future directions.

2. Preliminaries

2.1. The Ball–Plate Problem

The ball–plate problem is modeled by a mechanical system, and it is written as [28]

$$\left\{\begin{array}{c}\dot{x}=z(k-y)\hfill \\ \dot{y}=z(l+x)\hfill \\ \dot{z}=-kx-ly\hfill \end{array}\right.,k,l\in \mathbb{R}.$$

(1)

The mechanical system (1) has a Hamilton–Poisson structure characterized by the constants of motion given by

$$\begin{array}{c}H(x,y,z)={\displaystyle \frac{1}{2}}{(l+x)}^2+{\displaystyle \frac{1}{2}}{(k-y)}^2\hfill \\ C(x,y,z)={\displaystyle \frac{1}{2}}(x^2+y^2+z^2)\hfill \end{array},$$

(2)

with $H(x,y,z)$ the Hamiltonian and $C(x,y,z)$ the Casimir of the system.

The exact solution of Equations (1) and (4) is the intersection between the surfaces

$$\left\{\begin{array}{c}{(l+x)}^2+{(k-y)}^2={(l+x_0)}^2+{(k-y_0)}^2\hfill \\ x^2+y^2+z^2=x_0^2+y_0^2+z_0^2.\hfill \end{array}\right.,$$

(3)

with the initial conditions

$$x\left(0\right)=x_0,y\left(0\right)=y_0,z\left(0\right)=z_0.$$

(4)

2.2. Closed-Form Solutions

Three cases can be distinguished:

(a) F irstly, the physical parameters verify $k>0$, $l>0$.

By transformation,

$$\left\{\begin{array}{c}x=-l+R{\displaystyle \frac{1-w^2}{1+w^2}}\hfill \\ y=k-R{\displaystyle \frac{2w}{1+w^2}}\hfill \\ z=-2{\displaystyle \frac{\dot{w}}{1+w^2}}\hfill \end{array}\right.,$$

(5)

with $R=\sqrt{{(l+x_0)}^2+{(k-y_0)}^2}$, the closed-form solutions of the Equations (1) and (4) are obtained using the Hamilton–Poisson structure from Equation (2).

From the nonlinear initial value problem,

$$\left\{\begin{array}{c}-\ddot{w}(1+w^2)+2{\left(\dot{w}\right)}^{2}w+{\displaystyle \frac{kR}{2}}(1-w^4)-lRw(1+w^2)=0\hfill \\ w\left(0\right)={\displaystyle \frac{|k-y_0|}{|R+x_0+l|}},\dot{w}\left(0\right)=-{\displaystyle \frac{1}{2}}{z}_0(1+w^2\left(0\right)),\hfill \end{array}\right.$$

(6)

is obtained the unknown smooth function w, from Equation (5).

The set of the semi-analytical solutions of the Equation (6) will be denoted by

$${\mathcal{SA}}_1=\left\{\overline{w}\right|\overline{w}\text{semi-analytical} \mathrm{solution} \mathrm{of} \mathrm{the} \mathrm{Equation} \left(6\right)\}.$$

The substitution of $\overline{w}\in {\mathcal{SA}}_1$ into Equation (5) leads to semi-analytical closed-form solutions $\overline{x}$, $\overline{y}$, $\overline{z}$ of the Equations (1) and (4).

The set of the semi-analytical closed-form solutions of the Equations (1) and (4) will be denoted by

$$\begin{array}{c}{\mathcal{CFS}}_1=\{(\overline{x},\overline{y},\overline{z})|\overline{x},\overline{y},\overline{z}\text{semi-analytical} \text{closed-form} \mathrm{solution} \mathrm{of} \mathrm{the} \mathrm{Equations} \left(1\right) \mathrm{and} \left(4\right)\hfill \\ \mathrm{given} \mathrm{by} \mathrm{Equation}\left(5\right)\}.\hfill \end{array}$$

(b) In the second case the physical parameters satisfy $k<0$, $l>0$.

The closed-form solutions of the Equations (1) and (4) could be obtained using the following transformation:

$$\left\{\begin{array}{c}x=-l+R{\displaystyle \frac{1-w^2}{1+w^2}}\hfill \\ y=k+R{\displaystyle \frac{2w}{1+w^2}}\hfill \\ z=2{\displaystyle \frac{\dot{w}}{1+w^2}}\hfill \end{array}\right..$$

(7)

The unknown smooth function w from Equation (7) is solution of the nonlinear initial value problem:

$$\left\{\begin{array}{c}\ddot{w}(1+w^2)-2{\left(\dot{w}\right)}^{2}w+{\displaystyle \frac{kR}{2}}(1-w^4)+lRw(1+w^2)=0\hfill \\ w\left(0\right)={\displaystyle \frac{|k-y_0|}{|R+x_0+l|}},\dot{w}\left(0\right)={\displaystyle \frac{1}{2}}{z}_0(1+w^2\left(0\right)),\hfill \end{array}\right.$$

(8)

obtained from the third Equation (1).

The set of the semi-analytical solutions of the Equation (8) will be denoted by

$${\mathcal{SA}}_2=\left\{\overline{w}\right|\overline{w}\text{semi-analytical} \mathrm{solution} \mathrm{of} \mathrm{the} \mathrm{Equation} \left(8\right)\}.$$

The semi-analytical closed-form solutions $\overline{x}$, $\overline{y}$, $\overline{z}$ of the Equations (1) and (4) are obtained by substitution of the $\overline{w}\in {\mathcal{SA}}_2$ into Equation (7).

The set of the semi-analytical closed-form solutions of the Equations (1) and (4) will be denoted by

$$\begin{array}{c}{\mathcal{CFS}}_2=\{(\overline{x},\overline{y},\overline{z})|\overline{x},\overline{y},\overline{z}\text{semi-analytical} \text{closed-form} \mathrm{solution} \mathrm{of} \mathrm{the} \mathrm{Equations} \left(1\right) \mathrm{and} \left(4\right)\hfill \\ \mathrm{given} \mathrm{by} \mathrm{Equation}\left(7\right)\}.\hfill \end{array}$$

(c) The last case corresponds to the physical parameters $k\ne 0$, $l<0$.

By means of the transformation:

$$\left\{\begin{array}{c}x=-l+Rcos\left(w\right)\hfill \\ y=k+Rsin\left(w\right)\hfill \\ z=\dot{w}\hfill \end{array}\right.,$$

(9)

the closed-form solutions of the Equations (1) and (4) are built. The unknown smooth function w from Equation (9) is solution of the nonlinear initial value problem:

$$\left\{\begin{array}{c}\ddot{w}+kRcos\left(w\right)+lRsin\left(w\right)=0\hfill \\ w\left(0\right)=arctan\left({\displaystyle \frac{-k+y_0}{x_0+l}}\right),\dot{w}\left(0\right)=z_0,\hfill \end{array}\right.$$

(10)

obtained from the third Equation (1).

The set of the semi-analytical solutions of the Equation (10) will be denoted by

$${\mathcal{SA}}_3=\left\{\overline{w}\right|\overline{w}\text{semi-analytical} \mathrm{solution} \mathrm{of} \mathrm{the} \mathrm{Equation} \left(10\right)\}.$$

If $\overline{w}\in {\mathcal{SA}}_3$, then by substitution of $\overline{w}$ into Equation (9) are obtained the semi-analytical closed-form solutions $\overline{x}$, $\overline{y}$, $\overline{z}$ of the Equations (1) and (4).

The set of the semi-analytical closed-form solutions of the Equations (1) and (4) will be denoted by

$$\begin{array}{c}{\mathcal{CFS}}_3=\{(\overline{x},\overline{y},\overline{z})|\overline{x},\overline{y},\overline{z}\text{semi-analytical} \text{closed-form} \mathrm{solutions} \mathrm{of} \mathrm{the} \mathrm{Equations} \left(1\right) \mathrm{and} \left(4\right)\hfill \\ \mathrm{given} \mathrm{by} \mathrm{Equation}\left(9\right)\}.\hfill \end{array}$$

In the next section are presented some semi-analytical solutions of the Equations (6), (8) and (10) using the OPIM procedure and are denoted by $\overline{w}$ or ${\overline{w}}_{OPIM}$.

3. The Optimal Parametric Iteration Method

3.1. Preliminary

Let be the following second-order nonlinear differential equation:

$$\begin{array}{c}{\displaystyle \mathcal{L}\left[u\left(t\right)\right]+\mathcal{N}[t,u\left(t\right),\dot{u}\left(t\right),\ddot{u}\left(t\right)]-g\left(t\right)=0,t\in \mathcal{I}\subset \mathbb{R},}\hfill \end{array}$$

(11)

subject to the initial conditions

$$\begin{array}{c}{\displaystyle \mathcal{B}[u\left(0\right),\dot{u}\left(0\right)]=0,}\hfill \end{array}$$

(12)

with $\mathcal{L}$ a linear operator, $\mathcal{N}$ a nonlinear operator, $\mathcal{B}$ a boundary operator, g a known function, u an unknown smooth function depending on the independent variable t and $\dot{u}\left(t\right)={\displaystyle \frac{du}{dt}}$.

If an analytic function F is expanded in Taylor series it is obtained:

$$\begin{array}{c}F(t,u+\alpha ,\dot{u}+\beta ,\ddot{u}+\gamma )=F(t,u,\dot{u},\ddot{u})+{\displaystyle \frac{\alpha }{1!}}{F}_u(t,u,\dot{u},\ddot{u})+\hfill \\ +{\displaystyle \frac{\beta }{1!}}{F}_{\dot{u}}(t,u,\dot{u},\ddot{u})+{\displaystyle \frac{\gamma }{1!}}{F}_{\ddot{u}}(t,u,\dot{u},\ddot{u})+\dots ,\hfill \end{array}$$

(13)

where $F_u={\displaystyle \frac{\partial F}{\partial u}}$ denotes the partial derivative of the function F with respect to u, for $\alpha$, $\beta$, $\gamma$ real values. Instead of solving the nonlinear differential Equation (11), one can solve another equation, using Equation (40) and Optimal Parametric Iteration Method (OPIM), developed by Marinca et al. [20]:

$$\begin{array}{c}{\displaystyle \mathcal{L}\left[u_{n+1}\left(t\right)\right]+\mathcal{N}[t,u_n,{\dot{u}}_n,{\ddot{u}}_n]+{\alpha }_n(t,C_i){\mathcal{N}}_u[t,u_n,{\dot{u}}_n,{\ddot{u}}_n]+}\hfill \\ +{\beta }_n(t,C_j){\mathcal{N}}_{\dot{u}}[t,u_n,{\dot{u}}_n,{\ddot{u}}_n]+{\gamma }_n(t,C_k){\mathcal{N}}_{\ddot{u}}[t,u_n,{\dot{u}}_n,{\ddot{u}}_n]+\cdots -g\left(t\right)=0\hfill \\ {\displaystyle \mathcal{B}[u_{n+1}\left(t\right),{\dot{u}}_{n+1}\left(t\right)]=0,\mathrm{for}t=0}\hfill \end{array},n\ge 0,$$

(14)

where ${\alpha }_n(t,C_i)$, ${\beta }_n(t,C_j)$, and ${\gamma }_n(t,C_k)$ are auxiliary continuous functions; ${\mathcal{N}}_u={\displaystyle \frac{\partial \mathcal{N}}{\partial u}}$, ${\mathcal{N}}_{\dot{u}}={\displaystyle \frac{\partial \mathcal{N}}{\partial \dot{u}}}$, ${\mathcal{N}}_{\ddot{u}}={\displaystyle \frac{\partial \mathcal{N}}{\partial \ddot{u}}}$ (obtained from Taylor series expansion of the nonlinear operator $\mathcal{N}[t,u\left(t\right),\dot{u}\left(t\right),\ddot{u}\left(t\right)]$); $u_{n+1}\left(t\right)$ is the (n + 1)-th-order approximate solution of Equations (11) and (12), denoted by $\overline{u}\left(t\right)$ or ${\overline{u}}_{OPIM}\left(t\right)$ and $u_0\left(t\right)$ is the initial approximation, a solution of the linear differential problem:

$$\begin{array}{c}{\displaystyle \mathcal{L}\left[u_0\left(t\right)\right]-g\left(t\right)=0}\hfill \\ {\displaystyle \mathcal{B}[u_0\left(0\right),{\dot{u}}_0\left(0\right)]=0.}\hfill \end{array}$$

(15)

The unknown convergence-control parameters $C_i$, $C_j$, and $C_k$ that appears in Equation (14) will be optimally computed.

The $(n+1)$-order approximate solution of Equations (11) and (12) is well determined if the convergence-control parameters are known.

In OPIM, the linear operator $\mathcal{L}$ is arbitrarily chosen, not the physical parameters. There are situations when the choice of physical parameters give rise to chaotic behavior of the dynamic system. This happened in the case of choosing higher values of damping factor or for exceeding the optimal resonance conditions. Likewise in the case of the arbitrary choice of the initial conditions.

Let $u_0\left(t\right)$ be an initial approximation of Equation (15). The nonlinear operators $\mathcal{N}[t,u_0,{\dot{u}}_0,{\ddot{u}}_0]$, ${\mathcal{N}}_u[t,u_0,{\dot{u}}_0,{\ddot{u}}_0]$, ${\mathcal{N}}_{\dot{u}}[t,u_0,{\dot{u}}_0,{\ddot{u}}_0]$ and ${\mathcal{N}}_{\ddot{u}}[t,u_0,{\dot{u}}_0,{\ddot{u}}_0]$ that appear in Equation (14) have the form

$$\begin{array}{c}\sum _{i=1}^{n_{max}}{h}_i\left(t\right)g_i\left(t\right),\hfill \end{array}$$

(16)

where $n_{max}$ is a positive integer, and $h_i\left(t\right)$ and $g_i\left(t\right)$ are known functions that depend on $u_0\left(t\right)$.

For $u_{n+1}$ an $(n+1)$-order approximate solution of Equations (11) and (12), the validation of this procedure is highlighted by computing the residual function given by

$$\begin{array}{c}{\displaystyle \mathcal{R}\left(t\right)=\mathcal{L}\left[u_{n+1}\left(t\right)\right]+\mathcal{N}[t,u_{n+1}\left(t\right),{\dot{u}}_{n+1}\left(t\right),{\ddot{u}}_{n+1}\left(t\right)]-g\left(t\right),t\in \mathcal{I}\subset \mathbb{R},}\hfill \end{array}$$

(17)

such that ${\displaystyle \mathcal{R}\left(t\right)<<1}$, for all $t\in \mathcal{I}$.

Some mathematical notions as: OPIM sequence of the Equation (11), OPIM functions of Equation (11), $\epsilon$-approximate OPIM solutions of Equation (11), weak $\epsilon$-approximate OPIM solutions of Equation (11) on the real interval $(0,\infty )$ are defined in [29]. The existence of weak $\epsilon$-approximate OPIM solutions is built by the theorem presented in [29].

Remark 1. 

The integration of the Equation (14) produces secular terms of the form: $tcos\left({\omega }_{0}t\right)$, $tsin\left({\omega }_{0}t\right)$, $t^{2}cos\left({\omega }_{0}t\right)$, $t^{2}sin\left({\omega }_{0}t\right)$, $tcos\left(2{\omega }_{0}t\right)$, $tsin\left(2{\omega }_{0}t\right)$, and so on. For the nonlinear oscillator the secular terms that appear through integration generate the resonance phenomenon. Consequently, the secular terms have to be avoided.

3.2. Semi-Analytical Solutions via OPIM Technique

The applicability of the OPIM scheme for the Equation (6) using only one iteration is presented in details below.

$$\begin{array}{c}{\displaystyle \mathcal{L}\left[w\left(t\right)\right]=\ddot{w}+{\omega }_0^{2}w}\hfill \\ {\displaystyle \mathcal{N}[t,w\left(t\right),\dot{w}\left(t\right),\ddot{w}\left(t\right)]=-{\omega }_0^{2}w+\ddot{w}{w}^2-2{\left(\dot{w}\right)}^{2}w-\frac{kR}{2}(1-w^4)+lRu(1+w^2),t>0.}\hfill \end{array}$$

(18)

By means of the linear operator given by Equation (18), the initial approximation $w_0\left(t\right)$, solution of Equation (15) is

$$\begin{array}{c}{\displaystyle w_0\left(t\right)=B_{1}cos\left({\omega }_{0}t\right)+C_{1}sin\left({\omega }_{0}t\right),}\hfill \end{array}$$

(19)

with ${\omega }_0>0$, $B_1$, $C_1$ unknown parameters at this moment.

Using Equation (18), a simple computation yields the following expressions:

$$\begin{array}{c}{\displaystyle {\mathcal{N}}_w[t,w,\dot{w},\ddot{w}]=-{\omega }_0^2+2\ddot{w}w-2{\left(\dot{w}\right)}^2+2kR{w}^3+lR(1+w^2)+2lR{w}^2,}\hfill \\ {\mathcal{N}}_{\dot{w}}[t,w,\dot{w},\ddot{w}]=-4w\dot{w},{\mathcal{N}}_{\ddot{w}}[t,w,\dot{w},\ddot{w}]=w^2.\hfill \end{array}$$

(20)

For the Equation (14), there are more possibilities to choose the following auxiliary functions:

$$\begin{array}{c}{\displaystyle {\alpha }_n\left(t\right)=-\frac{1}{4}{B}_{1}cos\left({\omega }_{0}t\right)-\frac{1}{4}{C}_{1}sin\left({\omega }_{0}t\right),}\hfill \\ {\displaystyle {\beta }_n\left(t\right)=a_{1}cos\left({\omega }_{0}t\right)+b_{1}sin\left({\omega }_{0}t\right),}\hfill \\ {\displaystyle {\gamma }_n\left(t\right)=\sum _{j=3}^{N_{max}}{\tilde{B}}_{j}cos\left((2j+1){\omega }_{0}t\right)}\hfill \end{array}$$

(21)

or $\begin{array}{c}{\displaystyle {\alpha }_n\left(t\right)={\tilde{B}}_{2}cos\left(2{\omega }_{0}t\right)+{\tilde{C}}_{2}sin\left(2{\omega }_{0}t\right)}\hfill \end{array}$, ${\displaystyle {\beta }_n\left(t\right)={\tilde{B}}_{3}cos\left(3{\omega }_{0}t\right)+{\tilde{C}}_{3}sin\left(3{\omega }_{0}t\right)}$, ${\displaystyle {\gamma }_n\left(t\right)={\tilde{B}}_{4}cos\left(5{\omega }_{0}t\right)+{\tilde{C}}_{4}sin\left(5{\omega }_{0}t\right)}$, and so on.

By means of the Equations (19), (20) and (21) respectively, the following expression which appears in Equation (14) becomes:

$$\begin{array}{c}\mathcal{N}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\alpha }_n\left(t\right){\mathcal{N}}_w[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\beta }_n\left(t\right){\mathcal{N}}_{\dot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]={\displaystyle \frac{kR}{2}}+\left({\displaystyle \frac{3}{4}}{B}_{1}lR+\right.\hfill \\ +{\displaystyle \frac{3}{16}}{B}_1^{3}lR+{\displaystyle \frac{3}{16}}{B}_1{C}_1^{2}lR+B_1^2{b}_1{\omega }_0-2a_1{B}_1{C}_1{\omega }_0-b_1{C}_1^2{\omega }_0-{\displaystyle \frac{3}{4}}{B}_1{\omega }_0^2-\hfill \\ \left.-{\displaystyle \frac{3}{4}}{B}_1^3{\omega }_0^2-{\displaystyle \frac{3}{4}}{B}_1{C}_1^2{\omega }_0^2\right)cos\left({\omega }_{0}t\right)+\left({\displaystyle \frac{1}{16}}{B}_1^{3}lR-{\displaystyle \frac{3}{16}}{B}_1{C}_1^{2}lR-B_1^2{b}_1{\omega }_0-2a_1{B}_1{C}_1{\omega }_0+\right.\hfill \\ \left.+b_1{C}_1^2{\omega }_0+{\displaystyle \frac{1}{4}}{B}_1^3{\omega }_0^2-{\displaystyle \frac{3}{4}}{B}_1{C}_1^2{\omega }_0^2\right)cos\left(3{\omega }_{0}t\right)+\left({\displaystyle \frac{3}{4}}{C}_{1}lR+{\displaystyle \frac{3}{16}}{B}_1^2{C}_{1}lR+{\displaystyle \frac{3}{16}}{C}_1^{3}lR+\right.\hfill \\ \left.+a_1{B}_1^2{\omega }_0+2B_1{b}_1{C}_1{\omega }_0-a_1{C}_1^2{\omega }_0-{\displaystyle \frac{3}{4}}{C}_1{\omega }_0^2-{\displaystyle \frac{3}{4}}{B}_1^2{C}_1{\omega }_0^2-{\displaystyle \frac{3}{4}}{C}_1^3{\omega }_0^2\right)sin\left({\omega }_{0}t\right)+\hfill \\ +\left({\displaystyle \frac{3}{16}}{B}_1^2{C}_{1}lR-{\displaystyle \frac{1}{16}}{C}_1^{3}lR+a_1{B}_1^2{\omega }_0-2B_1{b}_1{C}_1{\omega }_0-a_1{C}_1^2{\omega }_0+\right.\hfill \\ \left.+{\displaystyle \frac{3}{4}}{B}_1^2{C}_1{\omega }_0^2-{\displaystyle \frac{1}{4}}{C}_1^3{\omega }_0^2\right)sin\left(3{\omega }_{0}t\right).\hfill \end{array}$$

(22)

This expression contains some terms depending on the elementary functions $cos\left({\omega }_{0}t\right)$ and $sin\left({\omega }_{0}t\right)$. By integration of the Equation (17), for $n=0$, the first-order approximate solution $w_1$ will contain secular terms of the form $tcos\left({\omega }_{0}t\right)$ and $tsin\left({\omega }_{0}t\right)$ which have to go to zero. Therefore, avoiding these secular terms in the Equation (22) involves:

$$\begin{array}{c}{a}_1=-{\displaystyle \frac{3C_1}{16(B_1^2+C_1^2){\omega }_0}}(-4lR-B_1^{2}lR-C_1^{2}lR+4{\omega }_0^2+4B_1^2{\omega }_0^2+4C_1^2{\omega }_0^2)\hfill \\ b_1={\displaystyle \frac{3}{16(B_1^2+C_1^2){\omega }_0}}(-4B_{1}lR-B_1^{3}lR-B_1{C}_1^{2}lR+4B_1{\omega }_0^2+4B_1^3{\omega }_0^2+4B_1{C}_1^2{\omega }_0^2)\hfill \end{array}.$$

(23)

For $N_{max}=3$, for example

$$\begin{array}{c}{\gamma }_n\left(t\right){\mathcal{N}}_{\ddot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]={\displaystyle \frac{1}{4}}{\tilde{B}}_3\left(B_1^2-C_1^2\right)cos\left(5{\omega }_{0}t\right)+{\displaystyle \frac{1}{4}}{\tilde{B}}_3\left(2B_1^2+2C_1^2\right)cos\left(7{\omega }_{0}t\right)+\hfill \\ +{\displaystyle \frac{1}{4}}{\tilde{B}}_3\left(B_1^2-C_1^2\right)cos\left(9{\omega }_{0}t\right)-{\displaystyle \frac{1}{2}}{\tilde{B}}_3{B}_1{C}_{1}sin\left(5{\omega }_{0}t\right)+{\displaystyle \frac{1}{2}}{\tilde{B}}_3{B}_1{C}_{1}sin\left(9{\omega }_{0}t\right).\hfill \end{array}$$

(24)

Thus, the expression $\mathcal{N}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\alpha }_n(t,C_i){\mathcal{N}}_w[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]++{\beta }_n(t,C_j){\mathcal{N}}_{\dot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\gamma }_n(t,C_k){\mathcal{N}}_{\ddot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]$ contains a linear combination of the elementary functions of the functions set

$$\left\{1,cos\left((2j+1){\omega }_{0}t\right),sin\left((2j+1){\omega }_{0}t\right)|j=\overline{1,4}\right\}.$$

In this way, by integration of the Equation (17), for $n=0$, the first-order approximate solution $w_1$ has the following form:

$$\begin{array}{c}{\displaystyle w_1\left(t\right)=D_0+\sum _{j=0}^4{D}_{j+1}cos\left((2j+1){\omega }_{0}t\right)+E_{j+1}sin\left((2j+1){\omega }_{0}t\right),}\hfill \end{array}$$

(25)

where the unknown control-convergence parameters $D_0$, $D_j$, $D_j$, $j=\overline{0,4}$ depend on the ${\omega }_0$, $B_1$, $C_1$, ${\tilde{B}}_3$ and will be optimally computed.

In the case when $N_{max}=4$ the expression ${\gamma }_n\left(t\right){\mathcal{N}}_{\ddot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]$ becomes:

$$\begin{array}{c}{\gamma }_n\left(t\right){\mathcal{N}}_{\ddot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]={\displaystyle \frac{1}{4}}{\tilde{B}}_3\left(B_1^2-C_1^2\right)cos\left(5{\omega }_{0}t\right)+{\displaystyle \frac{1}{4}}\left(2{\tilde{B}}_3{B}_1^2+2{\tilde{B}}_3{C}_1^2+{\tilde{B}}_4{B}_1^2-\right.\hfill \\ \left.-{\tilde{B}}_4{C}_1^2\right)cos\left(7{\omega }_{0}t\right)+{\displaystyle \frac{1}{4}}\left({\tilde{B}}_3{B}_1^2-{\tilde{B}}_3{C}_1^2+2{\tilde{B}}_4{B}_1^2+2{\tilde{B}}_4{C}_1^2\right)cos\left(9{\omega }_{0}t\right)+\hfill \\ +{\displaystyle \frac{1}{4}}{\tilde{B}}_4\left(B_1^2-C_1^2\right)cos\left(11{\omega }_{0}t\right)-{\displaystyle \frac{1}{2}}{\tilde{B}}_3{B}_1{C}_{1}sin\left(5{\omega }_{0}t\right)-{\displaystyle \frac{1}{2}}{\tilde{B}}_4{B}_1{C}_{1}sin\left(7{\omega }_{0}t\right)+\hfill \\ +{\displaystyle \frac{1}{2}}{\tilde{B}}_3{B}_1{C}_{1}sin\left(9{\omega }_{0}t\right)+{\displaystyle \frac{1}{2}}{\tilde{B}}_4{B}_1{C}_{1}sin\left(11{\omega }_{0}t\right).\hfill \end{array}$$

(26)

Thus, the expression $\mathcal{N}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\alpha }_n(t,C_i){\mathcal{N}}_w[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]++{\beta }_n(t,C_j){\mathcal{N}}_{\dot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\gamma }_n(t,C_k){\mathcal{N}}_{\ddot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]$ contains a linear combination of the elementary functions of the functions set

$$\left\{1,cos\left((2j+1){\omega }_{0}t\right),sin\left((2j+1){\omega }_{0}t\right)|j=\overline{1,5}\right\}.$$

As in the previous case, the first-order approximate solution $w_1$ is a linear combination of the elementary functions of the functions set

$$\left\{1,cos\left((2j+1){\omega }_{0}t\right),sin\left((2j+1){\omega }_{0}t\right)|j=\overline{0,5}\right\},$$

and so on.

Generally, for $N_{max}\in {\mathbb{N}}^{*}$, $N_{max}>4$ a fixed number the first approximation $w_1\left(t\right)$ is a linear combination of the functions set

$$\{1,cos\left((2j+1){\omega }_{0}t\right),sin\left((2j+1){\omega }_{0}t\right),j=0,1,\dots ,N_{max}+1\},$$

$$\begin{array}{c}{\displaystyle w_1\left(t\right)=D_0+\sum _{j=0}^{N_{max}+1}{D}_{j+1}cos\left((2j+1){\omega }_{0}t\right)+E_{j+1}sin\left((2j+1){\omega }_{0}t\right),}\hfill \end{array}$$

(27)

where the unknown control-convergence parameters $D_0$, $D_j$, $j=\overline{0,N_{max}+1}$ depend on the ${\omega }_0$, $B_1$, $C_1$ and will be optimally computed.

Thus, by using only one iteration, the OPIM solution is well determined as ${\overline{w}}_{OPIM}\left(t\right)=w_1\left(t\right)$ by Equation (27).

Furthermore, using two iterations, the OPIM solution is computed as ${\overline{w}}_{OPIM}\left(t\right)=w_2\left(t\right)$ by Equation (27), and so on.

For the first-order approximate solution $\overline{w}=w_1$ given by Equation (27) imposing the initial data from Equation (6) yields:

$$\begin{array}{c}{\displaystyle D_0=w\left(0\right)-\sum _{j=0}^{N_{max}+1}{D}_{j+1},{\omega }_0=\frac{\dot{w}\left(0\right)}{{\displaystyle \sum _{j=0}^{N_{max}+1}(2j+1)E_{j+1}}}}\hfill \end{array}.$$

(28)

In the second case of the Equation (8) by the same manner, the OPIM procedure is applied using only one iteration.

The linear and nonlinear operators, respectively, are:

$$\begin{array}{c}{\displaystyle \mathcal{L}\left[w\left(t\right)\right]=\ddot{w}\left(t\right)+{\omega }_0^{2}w\left(t\right)}\hfill \\ {\displaystyle \mathcal{N}[t,w\left(t\right),\dot{w}\left(t\right),\ddot{w}\left(t\right)]=-{\omega }_0^{2}w\left(t\right)+\ddot{w}{w}^2-2{\left(\dot{w}\right)}^{2}w+\frac{kR}{2}(1-w^4)+lRw(1+w^2),t>0.}\hfill \end{array}$$

(29)

Taking into consideration the linear operator given by Equation (29), the initial approximation $u_0\left(t\right)$, solution of Equation (15) is

$$\begin{array}{c}{\displaystyle w_0\left(t\right)=B_{1}cos\left({\omega }_{0}t\right)+C_{1}sin\left({\omega }_{0}t\right),}\hfill \end{array}$$

(30)

with ${\omega }_0>0$, $B_1$, $C_1$ unknown parameters at this moment.

From Equation (29), by a trivial computation, results the expressions:

$$\begin{array}{c}{\displaystyle {\mathcal{N}}_w[t,w,\dot{w},\ddot{w}]=-{\omega }_0^2+2\ddot{w}w-2{\left(\dot{w}\right)}^2-2kR{w}^3+lR(1+w^2)+2lR{w}^2,}\hfill \\ {\mathcal{N}}_{\dot{w}}[t,w,\dot{w},\ddot{w}]=-4w\dot{w},{\mathcal{N}}_{\ddot{w}}[t,w,\dot{w},\ddot{w}]=w^2.\hfill \end{array}$$

(31)

Returning to Equation (14), there are a lot of possibilities to choose the following auxiliary functions:

$$\begin{array}{c}{\displaystyle {\alpha }_n\left(t\right)=1,}\hfill \\ {\displaystyle {\beta }_n\left(t\right)=a_{1}cos\left({\omega }_{0}t\right)+b_{1}sin\left({\omega }_{0}t\right),}\hfill \\ {\displaystyle {\gamma }_n\left(t\right)=\sum _{j=3}^{N_{max}}{\tilde{B}}_{j}cos\left((2j+1){\omega }_{0}t\right)}\hfill \end{array}$$

(32)

or ${\displaystyle {\alpha }_n\left(t\right)={\tilde{B}}_{2}cos\left(2{\omega }_{0}t\right)}$, ${\displaystyle {\beta }_n\left(t\right)={\tilde{B}}_{3}cos\left(3{\omega }_{0}t\right)+{\tilde{C}}_{3}sin\left(3{\omega }_{0}t\right)}$, ${\displaystyle {\gamma }_n\left(t\right)={\tilde{B}}_{4}cos\left(5{\omega }_{0}t\right)+{\tilde{C}}_{4}cos\left(7{\omega }_{0}t\right)}$, and so on.

By means of the Equations (30), (31) and (32) respectively, the following expression from Equation (14) becomes:

$$\begin{array}{c}\mathcal{N}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\alpha }_n\left(t\right){\mathcal{N}}_w[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\beta }_n\left(t\right){\mathcal{N}}_{\dot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]=\hfill \\ ={\displaystyle \frac{kR}{2}}-{\displaystyle \frac{3}{16}}{B}_1^{4}kR-{\displaystyle \frac{3}{8}}{B}_1^2{C}_1^{2}kR-{\displaystyle \frac{3}{16}}{C}_1^{4}kR+{\displaystyle \frac{B_1^6{\omega }_0^4}{8}}+{\displaystyle \frac{3}{8}}{B}_1^4{C}_1^2{\omega }_0^4+{\displaystyle \frac{3}{8}}{B}_1^2{C}_1^4{\omega }_0^4+{\displaystyle \frac{C_1^6{\omega }_0^4}{8}}+\hfill \\ +lR+{\displaystyle \frac{3}{2}}{B}_1^{2}lR+{\displaystyle \frac{3}{2}}{C}_1^{2}lR-{\omega }_0^2+{\displaystyle \frac{1}{2}}{B}_1^4{\omega }_0^4+B_1^2{C}_1^2{\omega }_0^4+{\displaystyle \frac{1}{2}}{C}_1^4{\omega }_0^4+\left(B_{1}lR+{\displaystyle \frac{3}{4}}{B}_1^{3}lR+\right.\hfill \\ \left.+{\displaystyle \frac{3}{4}}{B}_1{C}_1^{2}lR+B_1^2{b}_1{\omega }_0-2a_1{B}_1{C}_1{\omega }_0-b_1{C}_1^2{\omega }_0-B_1{\omega }_0^2-{\displaystyle \frac{3}{2}}{B}_1^{3}kR-{\displaystyle \frac{3}{2}}{B}_1{C}_1^{2}kR\right)cos{\omega }_{0}t+\hfill \\ +\left(-{\displaystyle \frac{1}{4}}{B}_1^{4}k+{\displaystyle \frac{1}{4}}{C}_1^{4}kRR+{\displaystyle \frac{1}{16}}{B}_1^6{\omega }_0^4+{\displaystyle \frac{1}{16}}{B}_1^4{C}_1^2{\omega }_0^4-{\displaystyle \frac{1}{16}}{B}_1^2{C}_1^4{\omega }_0^4-{\displaystyle \frac{1}{16}}{C}_1^6{\omega }_0^4+{\displaystyle \frac{3}{2}}{B}_1^{2}lR-\right.\hfill \\ \left.-{\displaystyle \frac{3}{2}}{C}_1^{2}lR\right)cos2{\omega }_{0}t+\left({\displaystyle \frac{1}{4}}{B}_1^{3}lR-{\displaystyle \frac{3}{4}}{B}_1{C}_1^{2}lR-B_1^2{b}_1{\omega }_0-2a_1{B}_1{C}_1{\omega }_0+b_1{C}_1^2{\omega }_0-\right.\hfill \\ \left.-{\displaystyle \frac{1}{2}}{B}_1^{3}kR+{\displaystyle \frac{3}{2}}{B}_1{C}_1^{2}kR\right)cos3{\omega }_{0}t+\left(-{\displaystyle \frac{1}{16}}{B}_1^{4}kR+{\displaystyle \frac{3}{8}}{B}_1^2{C}_1^{2}kR-{\displaystyle \frac{1}{16}}{C}_1^{4}kR-{\displaystyle \frac{1}{8}}{B}_1^6{\omega }_0^4+\right.\hfill \\ \left.+{\displaystyle \frac{5}{8}}{B}_1^4{C}_1^2{\omega }_0^4+{\displaystyle \frac{5}{8}}{B}_1^2{C}_1^4{\omega }_0^4-{\displaystyle \frac{1}{8}}{C}_1^6{\omega }_0^4-{\displaystyle \frac{1}{2}}{B}_1^4{\omega }_0^4+3B_1^2{C}_1^2{\omega }_0^4-{\displaystyle \frac{1}{2}}{C}_1^4{\omega }_0^4\right)cos4{\omega }_{0}t+\hfill \\ +\left(-{\displaystyle \frac{1}{16}}{B}_1^6{\omega }_0^4+{\displaystyle \frac{15}{16}}{B}_1^4{C}_1^2{\omega }_0^4-{\displaystyle \frac{15}{16}}{B}_1^2{C}_1^4{\omega }_0^4+{\displaystyle \frac{1}{16}}{C}_1^6{\omega }_0^4\right)cos6{\omega }_{0}t+\left(C_{1}lR+{\displaystyle \frac{3}{4}}{B}_1^2{C}_{1}lR+\right.\hfill \\ \left.+{\displaystyle \frac{3}{4}}{C}_1^{3}lR+a_1{B}_1^2{\omega }_0+2B_1{b}_1{C}_1{\omega }_0-a_1{C}_1^2{\omega }_0-C_1{\omega }_0^2-{\displaystyle \frac{3}{2}}{B}_1^2{C}_{1}kR-{\displaystyle \frac{3}{2}}{C}_1^{3}kR\right)sin{\omega }_{0}t+\hfill \\ +\left(-{\displaystyle \frac{1}{2}}{B}_1^3{C}_{1}kR-{\displaystyle \frac{1}{2}}{B}_1{C}_1^{3}kR+{\displaystyle \frac{1}{8}}{B}_1^5{C}_1{\omega }_0^4+{\displaystyle \frac{1}{4}}{B}_1^3{C}_1^3{\omega }_0^4+{\displaystyle \frac{1}{8}}{B}_1{C}_1^5{\omega }_0^4+3B_1{C}_{1}lR\right)sin2{\omega }_{0}t+\hfill \\ +\left({\displaystyle \frac{3}{4}}{B}_1^2{C}_{1}lR-{\displaystyle \frac{1}{4}}{C}_1^{3}lR+a_1{B}_1^2{\omega }_0-2B_1{b}_1{C}_1{\omega }_0-a_1{C}_1^2{\omega }_0-{\displaystyle \frac{3}{2}}{B}_1^2{C}_{1}kR+\right.\hfill \\ \left.+{\displaystyle \frac{1}{2}}{C}_1^{3}kR\right)sin3{\omega }_{0}t+\left(-{\displaystyle \frac{1}{4}}{B}_1^3{C}_{1}kR+{\displaystyle \frac{1}{4}}{B}_1{C}_1^{3}kR-{\displaystyle \frac{1}{2}}{B}_1^5{C}_1{\omega }_0^4+{\displaystyle \frac{1}{2}}{B}_1{C}_1^5{\omega }_0^4-2B_1^3{C}_1{\omega }_0^4+\right.\hfill \\ \left.+2B_1{C}_1^3{\omega }_0^4\right)sin4{\omega }_{0}t+\left(-{\displaystyle \frac{3}{8}}{B}_1^5{C}_1{\omega }_0^4+{\displaystyle \frac{5}{4}}{B}_1^3{C}_1^3{\omega }_0^4-{\displaystyle \frac{3}{8}}{B}_1{C}_1^5{\omega }_0^4\right)sin6{\omega }_{0}t.\hfill \end{array}$$

(33)

This expression contains some terms depending on the elementary functions $cos\left({\omega }_{0}t\right)$ and $sin\left({\omega }_{0}t\right)$. By integration of the Equation (17), for $n=0$, the first-order approximate solution $w_1$ will contain secular terms of the form $tcos\left({\omega }_{0}t\right)$ and $tsin\left({\omega }_{0}t\right)$ which have to go to zero. Therefore, avoiding these terms in the Equation (33) involves:

$$\begin{array}{c}{a}_1={\displaystyle \frac{C_1}{4(B_1^2+C_1^2){\omega }_0}}(4lR+3B_1^{2}lR+3C_1^{2}lR-4{\omega }_0^2)\hfill \\ b_1=-{\displaystyle \frac{1}{4(B_1^2+C_1^2){\omega }_0}}(4B_{1}lR+3B_1^{3}lR+3B_1{C}_1^{2}lR-4B_1{\omega }_0^2)\hfill \end{array}.$$

(34)

For $N_{max}=3$

$$\begin{array}{c}{\gamma }_n\left(t\right){\mathcal{N}}_{\ddot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]={\displaystyle \frac{1}{4}}{\tilde{B}}_3\left(B_1^2-C_1^2\right)cos\left(5{\omega }_{0}t\right)+{\displaystyle \frac{1}{4}}{\tilde{B}}_3\left(2B_1^2+2C_1^2\right)cos\left(7{\omega }_{0}t\right)+\hfill \\ +{\displaystyle \frac{1}{4}}{\tilde{B}}_3\left(B_1^2-C_1^2\right)cos\left(9{\omega }_{0}t\right)-{\displaystyle \frac{1}{2}}{\tilde{B}}_3{B}_1{C}_{1}sin\left(5{\omega }_{0}t\right)+{\displaystyle \frac{1}{2}}{\tilde{B}}_3{B}_1{C}_{1}sin\left(9{\omega }_{0}t\right).\hfill \end{array}$$

(35)

Thus, the expression $\mathcal{N}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\alpha }_n(t,C_i){\mathcal{N}}_w[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]++{\beta }_n(t,C_j){\mathcal{N}}_{\dot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\gamma }_n(t,C_k){\mathcal{N}}_{\ddot{u}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]$ contain a linear combination of the elementary functions of the functions set

$$\left\{cos\left(i{\omega }_{0}t\right),sin\left(i{\omega }_{0}t\right)|i=\overline{1,6}\right\}\cup \left\{1,cos\left((2j+1){\omega }_{0}t\right),sin\left((2j+1){\omega }_{0}t\right)|j=\overline{3,4}\right\}.$$

In this way, by integration of the Equation (17), for $n=0$, the first-order approximate solution $w_1$ has the following form:

$$\begin{array}{c}{\displaystyle w_1\left(t\right)=D_0+\sum _{i=1}^6{D}_{i}cos\left(i\right){\omega }_{0}t)+E_{i}sin\left(i{\omega }_{0}t\right)+}\hfill \\ {\displaystyle +\sum _{j=3}^4{D}_{j+4}cos\left((2j+1){\omega }_{0}t\right)+E_{j+4}sin\left((2j+1){\omega }_{0}t\right),}\hfill \end{array}$$

(36)

with unknown control-convergence parameters $D_j$, $j=\overline{0,8}$ depending on the ${\omega }_0$, $B_1$, $C_1$, ${\tilde{B}}_3$ and will be optimally computed.

In the case when $N_{max}=4$ the expression ${\gamma }_n\left(t\right){\mathcal{N}}_{\ddot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]$ become:

$$\begin{array}{c}{\gamma }_n\left(t\right){\mathcal{N}}_{\ddot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]={\displaystyle \frac{1}{4}}{\tilde{B}}_3\left(B_1^2-C_1^2\right)cos\left(5{\omega }_{0}t\right)+{\displaystyle \frac{1}{4}}\left(2{\tilde{B}}_3{B}_1^2+2{\tilde{B}}_3{C}_1^2+{\tilde{B}}_4{B}_1^2-\right.\hfill \\ \left.-{\tilde{B}}_4{C}_1^2\right)cos\left(7{\omega }_{0}t\right)+{\displaystyle \frac{1}{4}}\left({\tilde{B}}_3{B}_1^2-{\tilde{B}}_3{C}_1^2+2{\tilde{B}}_4{B}_1^2+2{\tilde{B}}_4{C}_1^2\right)cos\left(9{\omega }_{0}t\right)+\hfill \\ +{\displaystyle \frac{1}{4}}{\tilde{B}}_4\left(B_1^2-C_1^2\right)cos\left(11{\omega }_{0}t\right)-{\displaystyle \frac{1}{2}}{\tilde{B}}_3{B}_1{C}_{1}sin\left(5{\omega }_{0}t\right)-{\displaystyle \frac{1}{2}}{\tilde{B}}_4{B}_1{C}_{1}sin\left(7{\omega }_{0}t\right)+\hfill \\ +{\displaystyle \frac{1}{2}}{\tilde{B}}_3{B}_1{C}_{1}sin\left(9{\omega }_{0}t\right)+{\displaystyle \frac{1}{2}}{\tilde{B}}_4{B}_1{C}_{1}sin\left(11{\omega }_{0}t\right).\hfill \end{array}$$

(37)

Thus, the expression $\mathcal{N}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\alpha }_n(t,C_i){\mathcal{N}}_w[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]++{\beta }_n(t,C_j){\mathcal{N}}_{\dot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\gamma }_n(t,C_k){\mathcal{N}}_{\ddot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]$ contain a linear combination of the elementary functions of the functions set

$$\left\{cos\left(i{\omega }_{0}t\right),sin\left(i{\omega }_{0}t\right)|i=\overline{1,6}\right\}\cup \left\{1,cos\left((2j+1){\omega }_{0}t\right),sin\left((2j+1){\omega }_{0}t\right)|j=\overline{3,5}\right\}.$$

As in the previous case, the first-order approximate solution $w_1$ is a linear combination of the elementary functions of the functions set

$$\left\{cos\left(i{\omega }_{0}t\right),sin\left(i{\omega }_{0}t\right)|i=\overline{1,6}\right\}\cup \left\{1,cos\left((2j+1){\omega }_{0}t\right),sin\left((2j+1){\omega }_{0}t\right)|j=\overline{3,5}\right\},$$

and so on.

Generally, for $N_{max}\in {\mathbb{N}}^{*}$, $N_{max}>4$ a fixed number, the first approximation $w_1\left(t\right)$ is a linear combination of the functions set

$$\left\{cos\left(i{\omega }_{0}t\right),sin\left(i{\omega }_{0}t\right)|i=\overline{1,6}\right\}\cup \left\{1,cos\left((2j+1){\omega }_{0}t\right),sin\left((2j+1){\omega }_{0}t\right)|j=\overline{3,N_{max}+1}\right\}$$

$$\begin{array}{c}{\displaystyle w_1\left(t\right)=D_0+\sum _{i=1}^6{D}_{i}cos\left(i\right){\omega }_{0}t)+E_{i}sin\left(i{\omega }_{0}t\right)+}\hfill \\ {\displaystyle +\sum _{j=3}^{N_{max}+1}{D}_{j+4}cos\left((2j+1){\omega }_{0}t\right)+E_{j+4}sin\left((2j+1){\omega }_{0}t\right),}\hfill \end{array}$$

(38)

where the unknown control-convergence parameters $D_0$, $D_j$, $D_j$, $j=\overline{1,N_{max}+1}$ depend on the ${\omega }_0$, $B_1$, $C_1$, ${\tilde{B}}_i$, $i=\overline{1,N_{max}}$ and will be optimally computed.

Thus, by using only one iteration, the OPIM solution is well determined as ${\overline{w}}_{OPIM}\left(t\right)=w_1\left(t\right)$ by Equation (38).

Furthermore, using two iterations, the OPIM solution is computed as ${\overline{w}}_{OPIM}\left(t\right)=w_2\left(t\right)$ by Equation (38), and so on.

Let $\overline{w}=w_1$ be a first-order approximate solution given by Equation (38). Imposing the initial data from Equation (8) yields:

$$\begin{array}{c}{\displaystyle D_0=w\left(0\right)-\sum _{j=1}^{N_{max}+5}{D}_j,{\omega }_0=\frac{\dot{w}\left(0\right)}{{\displaystyle \sum _{i=1}^{6}i{E}_i+\sum _{j=3}^{N_{max}+1}(2j+1)E_{j+4}}}}\hfill \end{array}.$$

(39)

The applicability of the OPIM procedure for the nonlinear differential problem given by Equation (10) using only one iteration is presented in details below.

For the third case $k\ne 0$, $l<0$ the following well-known series expansions are used:

$$\begin{array}{c}{\displaystyle cos\left(w\right)=\sum _{j=0}^{\infty }\frac{{(-1)}^j}{\left(2j\right)!}{w}^{2j}=1-\frac{1}{2!}{w}^2+\frac{1}{4!}{w}^4-\frac{1}{6!}{w}^6+\dots +\frac{{(-1)}^{N_{max}}}{\left(2N_{max}\right)!}{w}^{2N_{max}}+\dots ,}\hfill \\ {\displaystyle sin\left(w\right)=\sum _{j=0}^{\infty }\frac{{(-1)}^j}{(2j+1)!}{w}^{2j+1}=w-\frac{1}{3!}{w}^3+\frac{1}{5!}{w}^5-\frac{1}{7!}{w}^7+\dots +\frac{{(-1)}^{N_{max}}}{(2N_{max}+1)!}{w}^{2N_{max}+1}+\dots .}\hfill \end{array}.$$

(40)

The linear and nonlinear operators could be chosen as:

$$\begin{array}{c}{\displaystyle \mathcal{L}\left[w\left(t\right)\right]=\ddot{w}+{\omega }_0^{2}w\left(t\right)}\hfill \\ {\displaystyle \mathcal{N}[t,w\left(t\right),\dot{w}\left(t\right),\ddot{w}\left(t\right)]=-{\omega }_0^{2}w+kRcos\left(w\right)+lRsin\left(w\right),t>0.}\hfill \end{array}$$

(41)

Taking into consideration the linear operator given by Equation (41), the initial approximation $w_0\left(t\right)$, solution of Equation (15) is

$$\begin{array}{c}{\displaystyle w_0\left(t\right)=B_{1}cos\left({\omega }_{0}t\right)+C_{1}sin\left({\omega }_{0}t\right),}\hfill \end{array}$$

(42)

with ${\omega }_0>0$, $B_1$, $C_1$ unknown parameters at this moment.

Using Equation (41), a simple computation yields the following expressions:

$$\begin{array}{c}{\displaystyle {\mathcal{N}}_w[t,w,\dot{w},\ddot{w}]=-{\omega }_0^2-kRsin\left(w\right)+lRcos\left(w\right),}\hfill \\ {\mathcal{N}}_{\dot{w}}[t,w,\dot{w},\ddot{w}]=0,{\mathcal{N}}_{\ddot{w}}[t,w,\dot{w},\ddot{w}]=0.\hfill \end{array}$$

(43)

In this case the following auxiliary functions are proposed:

$$\begin{array}{c}{\displaystyle {\alpha }_n\left(t\right)=a_{1}cos\left({\omega }_{0}t\right)+b_1,}\hfill \\ {\displaystyle {\beta }_n\left(t\right)=1,}\hfill \\ {\displaystyle {\gamma }_n\left(t\right)=1}\hfill \end{array}$$

(44)

or ${\displaystyle {\alpha }_n\left(t\right)={\tilde{B}}_{2}cos\left(2{\omega }_{0}t\right)}$, ${\displaystyle {\beta }_n\left(t\right)={\tilde{B}}_{3}cos\left(3{\omega }_{0}t\right)+{\tilde{C}}_{3}sin\left(3{\omega }_{0}t\right)}$, ${\displaystyle {\gamma }_n\left(t\right)={\tilde{B}}_{4}cos\left(5{\omega }_{0}t\right)+{\tilde{C}}_{4}cos\left(7{\omega }_{0}t\right)}$, and so on.

A semi-analytical solution of the Equation (10) can be obtain by OPIM procedure choosing an arbitrary value of the index $N_{max}$ in the Equation (40).

Thus, for $N_{max}=2$, using Equations (40) and (43), the expression $\mathcal{N}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\alpha }_n\left(t\right){\mathcal{N}}_w[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]$ has the form:

$$\begin{array}{c}\mathcal{N}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\alpha }_n\left(t\right){\mathcal{N}}_w[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]=D_0+D_{1}cos\left({\omega }_{0}t\right)+E_{1}sin\left({\omega }_{0}t\right)+\hfill \\ +D_{2}cos\left(2{\omega }_{0}t\right)+E_{2}sin\left(2{\omega }_{0}t\right)+D_{3}cos\left(3{\omega }_{0}t\right)+E_{3}sin\left(3{\omega }_{0}t\right)+D_{4}cos\left(4{\omega }_{0}t\right)+\hfill \\ +E_{4}sin\left(4{\omega }_{0}t\right)+D_{5}cos\left(5{\omega }_{0}t\right)+E_{5}sin\left(5{\omega }_{0}t\right)+D_{6}cos\left(6{\omega }_{0}t\right)+E_{6}sin\left(6{\omega }_{0}t\right),\hfill \end{array}$$

(45)

where

$\begin{array}{c}{D}_0=kR-{\displaystyle \frac{1}{2}}{a}_1{B}_{1}kR-{\displaystyle \frac{1}{4}}{B}_1^{2}kR+{\displaystyle \frac{1}{16}}{a}_1{B}_1^{3}kR+{\displaystyle \frac{1}{64}}{B}_1^{4}kR-{\displaystyle \frac{1}{384}}{a}_1{B}_1^{5}kR-{\displaystyle \frac{1}{4}}{C}_1^{2}kR+\hfill \\ +{\displaystyle \frac{1}{16}}{a}_1{B}_1{C}_1^{2}kR+{\displaystyle \frac{1}{32}}{B}_1^2{C}_1^{2}kR-{\displaystyle \frac{1}{192}}{a}_1{B}_1^3{C}_1^{2}kR+{\displaystyle \frac{1}{64}}{C}_1^{4}kR-{\displaystyle \frac{1}{384}}{a}_1{B}_1{C}_1^{4}kR+b_{1}lR-\hfill \\ -{\displaystyle \frac{1}{4}}{B}_1^2{b}_{1}lR+{\displaystyle \frac{1}{64}}{B}_1^4{b}_{1}lR-{\displaystyle \frac{1}{4}}{b}_1{C}_1^{2}lR+{\displaystyle \frac{1}{32}}{B}_1^2{b}_1{C}_1^{2}lR+{\displaystyle \frac{1}{64}}{b}_1{C}_1^{4}lR-b_1{\omega }_0^2,\hfill \end{array}$

$\begin{array}{c}{D}_1=-B_1{b}_{1}kR+{\displaystyle \frac{1}{8}}{B}_1^3{b}_{1}kR-{\displaystyle \frac{1}{192}}{B}_1^5{b}_{1}kR+{\displaystyle \frac{1}{8}}{B}_1{b}_1{C}_1^{2}kR-{\displaystyle \frac{1}{96}}{B}_1^3{b}_1{C}_1^{2}kR-\hfill \\ -{\displaystyle \frac{1}{192}}{B}_1{b}_1{C}_1^{4}kR+a_{1}lR+B_{1}lR-{\displaystyle \frac{3}{8}}{a}_1{B}_1^{2}lR-{\displaystyle \frac{1}{8}}{B}_1^{3}lR+{\displaystyle \frac{5}{192}}{a}_1{B}_1^{4}lR+{\displaystyle \frac{1}{192}}{B}_1^{5}lR-{\displaystyle \frac{1}{8}}{a}_1{C}_1^{2}lR-\hfill \\ -{\displaystyle \frac{1}{8}}{B}_1{C}_1^{2}lR+{\displaystyle \frac{1}{32}}{a}_1{B}_1^2{C}_1^{2}lR+{\displaystyle \frac{1}{96}}{B}_1^3{C}_1^{2}lR+{\displaystyle \frac{1}{192}}{a}_1{C}_1^{4}lR+{\displaystyle \frac{1}{192}}{B}_1{C}_1^{4}lR-a_1{\omega }_0^2-B_1{\omega }_0^2,\hfill \end{array}$

$\begin{array}{c}{D}_2=-{\displaystyle \frac{1}{2}}{a}_1{B}_{1}kR-{\displaystyle \frac{1}{4}}{B}_1^{2}kR+{\displaystyle \frac{1}{12}}{a}_1{B}_1^{3}kR+{\displaystyle \frac{1}{48}}{B}_1^{4}kR-{\displaystyle \frac{1}{256}}{a}_1{B}_1^{5}kR+{\displaystyle \frac{1}{4}}{C}_1^{2}kR-\hfill \\ -{\displaystyle \frac{1}{384}}{a}_1{B}_1^3{C}_1^{2}kR-{\displaystyle \frac{1}{48}}{C}_1^{4}kR+{\displaystyle \frac{1}{768}}{a}_1{B}_1{C}_1^{4}kR-{\displaystyle \frac{1}{4}}{B}_1^2{b}_{1}lR+{\displaystyle \frac{1}{48}}{B}_1^4{b}_{1}lR+{\displaystyle \frac{1}{4}}{b}_1{C}_1^{2}lR-{\displaystyle \frac{1}{48}}{b}_1{C}_1^{4}lR,\hfill \end{array}$

$\begin{array}{c}{D}_3={\displaystyle \frac{1}{24}}{B}_1^3{b}_{1}kR-{\displaystyle \frac{1}{384}}{B}_1^5{b}_{1}kR-{\displaystyle \frac{1}{8}}{B}_1{b}_1{C}_1^{2}kR+{\displaystyle \frac{1}{192}}{B}_1^3{b}_1{C}_1^{2}kR+{\displaystyle \frac{1}{128}}{B}_1{b}_1{C}_1^{4}kR-\hfill \\ -{\displaystyle \frac{1}{8}}{a}_1{B}_1^{2}lR-{\displaystyle \frac{1}{24}}{B}_1^{3}lR+{\displaystyle \frac{5}{384}}{a}_1{B}_1^{4}lR+{\displaystyle \frac{1}{384}}{B}_1^{5}lR+{\displaystyle \frac{1}{8}}{a}_1{C}_1^{2}lR+{\displaystyle \frac{1}{8}}{B}_1{C}_1^{2}lR-{\displaystyle \frac{1}{64}}{a}_1{B}_1^2{C}_1^{2}lR-\hfill \\ -{\displaystyle \frac{1}{192}}{B}_1^3{C}_1^{2}lR-{\displaystyle \frac{1}{128}}{a}_1{C}_1^{4}lR-{\displaystyle \frac{1}{128}}{B}_1{C}_1^{4}lR,\hfill \end{array}$

$\begin{array}{c}{D}_4={\displaystyle \frac{1}{48}}{a}_1{B}_1^{3}kR+{\displaystyle \frac{1}{192}}{B}_1^{4}kR-{\displaystyle \frac{1}{640}}{a}_1{B}_1^{5}kR-{\displaystyle \frac{1}{16}}{a}_1{B}_1{C}_1^{2}kR-{\displaystyle \frac{1}{32}}{B}_1^2{C}_1^{2}kR+{\displaystyle \frac{1}{192}}{a}_1{B}_1^3{C}_1^{2}kR+\hfill \\ +{\displaystyle \frac{1}{192}}{C}_1^{4}kR+{\displaystyle \frac{1}{384}}{a}_1{B}_1{C}_1^{4}kR+{\displaystyle \frac{1}{192}}{B}_1^4{b}_{1}lR-{\displaystyle \frac{1}{32}}{B}_1^2{b}_1{C}_1^{2}lR+{\displaystyle \frac{1}{192}}{b}_1{C}_1^{4}lR,\hfill \end{array}$

$\begin{array}{c}{D}_5=-{\displaystyle \frac{1}{1920}}{B}_1^5{b}_{1}kR+{\displaystyle \frac{1}{192}}{B}_1^3{b}_1{C}_1^{2}kR-{\displaystyle \frac{1}{384}}{B}_1{b}_1{C}_1^{4}kR+{\displaystyle \frac{1}{384}}{a}_1{B}_1^{4}lR+{\displaystyle \frac{1}{1920}}{B}_1^{5}lR-\hfill \\ -{\displaystyle \frac{1}{64}}{a}_1{B}_1^2{C}_1^{2}lR-{\displaystyle \frac{1}{192}}{B}_1^3{C}_1^{2}lR+{\displaystyle \frac{1}{384}}{a}_1{C}_1^{4}lR+{\displaystyle \frac{1}{384}}{B}_1{C}_1^{4}lR,\hfill \end{array}$

$\begin{array}{c}{D}_6=-{\displaystyle \frac{1}{3840}}{a}_1{B}_1^{5}kR+{\displaystyle \frac{1}{384}}{a}_1{B}_1^3{C}_1^{2}kR-{\displaystyle \frac{1}{768}}{a}_1{B}_1{C}_1^{4}kR,\hfill \end{array}$

$\begin{array}{c}{E}_1=-b_1{C}_{1}kR+{\displaystyle \frac{1}{8}}{B}_1^2{b}_1{C}_{1}kR-{\displaystyle \frac{1}{192}}{B}_1^4{b}_1{C}_{1}kR+{\displaystyle \frac{1}{8}}{b}_1{C}_1^{3}kR-{\displaystyle \frac{1}{96}}{B}_1^2{b}_1{C}_1^{3}kR-\hfill \\ -{\displaystyle \frac{1}{192}}{b}_1{C}_1^{5}kR+C_{1}lR-{\displaystyle \frac{1}{4}}{a}_1{B}_1{C}_{1}lR-{\displaystyle \frac{1}{8}}{B}_1^2{C}_{1}lR+{\displaystyle \frac{1}{48}}{a}_1{B}_1^3{C}_{1}lR+{\displaystyle \frac{1}{192}}{B}_1^4{C}_{1}lR-\hfill \\ -{\displaystyle \frac{1}{8}}{C}_1^{3}lR+{\displaystyle \frac{1}{48}}{a}_1{B}_1{C}_1^{3}lR+{\displaystyle \frac{1}{96}}{B}_1^2{C}_1^{3}lR+{\displaystyle \frac{1}{192}}{C}_1^{5}lR-C_1{\omega }_0^2,\hfill \end{array}$

$\begin{array}{c}{E}_2=-{\displaystyle \frac{1}{2}}{a}_1{C}_{1}kR-{\displaystyle \frac{1}{2}}{B}_1{C}_{1}kR+{\displaystyle \frac{1}{8}}{a}_1{B}_1^2{C}_{1}kR+{\displaystyle \frac{1}{24}}{B}_1^3{C}_{1}kR-{\displaystyle \frac{5}{768}}{a}_1{B}_1^4{C}_{1}kR+{\displaystyle \frac{1}{24}}{a}_1{C}_1^{3}kR+\hfill \\ +{\displaystyle \frac{1}{24}}{B}_1{C}_1^{3}kR-{\displaystyle \frac{1}{128}}{a}_1{B}_1^2{C}_1^{3}kR-{\displaystyle \frac{1}{768}}{a}_1{C}_1^{5}kR-{\displaystyle \frac{1}{2}}{B}_1{b}_1{C}_{1}lR+{\displaystyle \frac{1}{24}}{B}_1^3{b}_1{C}_{1}lR+{\displaystyle \frac{1}{24}}{B}_1{b}_1{C}_1^{3}lR,\hfill \end{array}$

$\begin{array}{c}{E}_3={\displaystyle \frac{1}{8}}{B}_1^2{b}_1{C}_{1}kR-{\displaystyle \frac{1}{128}}{B}_1^4{b}_1{C}_{1}kR-{\displaystyle \frac{1}{24}}{b}_1{C}_1^{3}kR-{\displaystyle \frac{1}{192}}{B}_1^2{b}_1{C}_1^{3}kR+{\displaystyle \frac{1}{384}}{b}_1{C}_1^{5}kR-\hfill \\ -{\displaystyle \frac{1}{4}}{a}_1{B}_1{C}_{1}lR-{\displaystyle \frac{1}{8}}{B}_1^2{C}_{1}lR+{\displaystyle \frac{1}{32}}{a}_1{B}_1^3{C}_{1}lR+{\displaystyle \frac{1}{128}}{B}_1^4{C}_{1}lR+{\displaystyle \frac{1}{24}}{C}_1^{3}lR+{\displaystyle \frac{1}{96}}{a}_1{B}_1{C}_1^{3}lR+\hfill \\ +{\displaystyle \frac{1}{192}}{B}_1^2{C}_1^{3}lR-{\displaystyle \frac{1}{384}}{C}_1^{5}lR,\hfill \end{array}$

$\begin{array}{c}{E}_4={\displaystyle \frac{1}{16}}{a}_1{B}_1^2{C}_{1}kR+{\displaystyle \frac{1}{48}}{B}_1^3{C}_{1}kR-{\displaystyle \frac{1}{192}}{a}_1{B}_1^4{C}_{1}kR-{\displaystyle \frac{1}{48}}{a}_1{C}_1^{3}kR-{\displaystyle \frac{1}{48}}{B}_1{C}_1^{3}kR+{\displaystyle \frac{1}{960}}{a}_1{C}_1^{5}kR+\hfill \\ +{\displaystyle \frac{1}{48}}{B}_1^3{b}_1{C}_{1}lR-{\displaystyle \frac{1}{48}}{B}_1{b}_1{C}_1^{3}lR,\hfill \end{array}$

$\begin{array}{c}{E}_5=-{\displaystyle \frac{1}{384}}{B}_1^4{b}_1{C}_{1}kR+{\displaystyle \frac{1}{192}}{B}_1^2{b}_1{C}_1^{3}kR-{\displaystyle \frac{1}{1920}}{b}_1{C}_1^{5}kR+{\displaystyle \frac{1}{96}}{a}_1{B}_1^3{C}_{1}lR+{\displaystyle \frac{1}{384}}{B}_1^4{C}_{1}lR-\hfill \\ -{\displaystyle \frac{1}{96}}{a}_1{B}_1{C}_1^{3}lR-{\displaystyle \frac{1}{192}}{B}_1^2{C}_1^{3}lR+{\displaystyle \frac{1}{1920}}{C}_1^{5}lR,\hfill \end{array}$

$\begin{array}{c}{E}_6=-{\displaystyle \frac{1}{768}}{a}_1{B}_1^4{C}_{1}kR+{\displaystyle \frac{1}{384}}{a}_1{B}_1^2{C}_1^{3}kR-{\displaystyle \frac{1}{3840}}{a}_1{C}_1^{5}kR.\hfill \end{array}$

Imposing the conditions $D_1=0$ and $E_1=0$ (by avoiding the secular terms in $cos\left({\omega }_{0}t\right)$ and $sin\left({\omega }_{0}t\right)$ ) involves:

$$\begin{array}{c}{a}_1=0\hfill \\ b_1=-{\displaystyle \frac{-192lR+24B_1^{2}lR-B_1^{4}lR+24C_1^{2}lR-2B_1^2{C}_1^{2}lR-C_1^{4}lR+192{\omega }_0^2}{(192-24B_1^2+B_1^4-24C_1^2+2B_1^2{C}_1^2+C_1^4)kR}}\hfill \end{array}.$$

(46)

Therefore, the expression $\mathcal{N}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\alpha }_n\left(t\right){\mathcal{N}}_w[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\beta }_n\left(t\right){\mathcal{N}}_{\dot{w}}$$[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]+{\gamma }_n\left(t\right){\mathcal{N}}_{\ddot{w}}[t,w_0,{\dot{w}}_0,{\ddot{w}}_0]$ becomes a linear combination of the functions set

$$\left\{1,cos\left(j{\omega }_{0}t\right),sin\left(j{\omega }_{0}t\right)|j=\overline{2,6}\right\}.$$

Consequently, the first-order approximate solution $w_1$, solution of the Equation (14) becomes a linear combination of the functions set

$$\left\{1,cos\left(j{\omega }_{0}t\right),sin\left(j{\omega }_{0}t\right)|j=\overline{1,6}\right\}.$$

Analogously, for $N_{max}=3$ the first-order approximate solution $w_1$, solution of the Equation (14) becomes a linear combination of the functions set

$$\left\{1,cos\left(j{\omega }_{0}t\right),sin\left(j{\omega }_{0}t\right)|j=\overline{1,8}\right\},$$

and so on.

For $N_{max}>3$ a fixed arbitrary number the first-order approximate solution $w_1$, solution of the Equation (14) has the form

$$\begin{array}{c}{\displaystyle w_1\left(t\right)=D_0+\sum _{j=1}^{2N_{max}+2}{D}_{j}cos\left(j{\omega }_{0}t\right)+E_{j}sin\left(j{\omega }_{0}t\right),}\hfill \end{array}$$

(47)

where the unknown control-convergence parameters $D_0$, $D_j$, $E_j$, $j=\overline{1,2N_{max}+2}$ depend on the ${\omega }_0$, $B_1$, $C_1$ and will be optimally computed.

If $\overline{w}=w_1$ is a first-order approximate solution given by Equation (47), then imposing the initial data from Equation (10) yields:

$$\begin{array}{c}{\displaystyle D_0=w\left(0\right)-\sum _{j=1}^{2N_{max}+2}{D}_j,{\omega }_0=\frac{\dot{w}\left(0\right)}{{\displaystyle \sum _{j=1}^{2N_{max}+2}j{E}_j}}}\hfill \end{array}.$$

(48)

4. Numerical versus OPIM Results

In this section are presented in details the numerical versus OPIM solutions for the ball–plate problem. The studied problem admits periodic solutions independently of initial conditions as it was proved in [28].

In the considered cases the numerical solutions are obtained via the 4th-order Runge-Kutta method.

4.1. Case 1: $k>0$, $l>0$

To highlight the accuracy of the obtained analytical solutions in

Table 1. Values of the absolute errors: ${ϵ}_w=|w_{numerical}-{\overline{w}}_{OPIM}|$ for initial data: $x_0=0.45$, $y_0=0.25$, $z_0=0.75$, physical parameters: $k=0.85$, $l=0.5$ and four values of the index $N_{max}$; ${\overline{w}}_{OPIM}\in {\mathcal{SA}}_1$ using Equations (27), (28) and (A1)–(A3).

t	${\mathit{ϵ}}_{\mathit{w}}:$ ${\mathit{N}}_{\mathit{max}}\mathbf{=}\mathbf{13}$	${\mathit{ϵ}}_{\mathit{w}}:$ ${\mathit{N}}_{\mathit{max}}\mathbf{=}\mathbf{18}$	${\mathit{ϵ}}_{\mathit{w}}:$ ${\mathit{N}}_{\mathit{max}}\mathbf{=}\mathbf{23}$	${\mathit{ϵ}}_{\mathit{w}}:$ ${\mathit{N}}_{\mathit{max}}\mathbf{=}\mathbf{28}$
0	1.2351 $\times {10}^{-13}$	1.8091 $\times {10}^{-12}$	1.8596 $\times {10}^{-14}$	6.9194 $\times {10}^{-13}$
2/3	2.6619 $\times {10}^{-5}$	4.2138 $\times {10}^{-7}$	8.1428 $\times {10}^{-8}$	2.9801 $\times {10}^{-9}$
4/3	2.8862 $\times {10}^{-5}$	6.4069 $\times {10}^{-7}$	7.3662 $\times {10}^{-8}$	2.5618 $\times {10}^{-9}$
2	2.2258 $\times {10}^{-5}$	3.5944 $\times {10}^{-7}$	5.6622 $\times {10}^{-8}$	2.1801 $\times {10}^{-9}$
8/3	1.1788 $\times {10}^{-6}$	2.2812 $\times {10}^{-7}$	7.1642 $\times {10}^{-8}$	2.5417 $\times {10}^{-9}$
10/3	4.5917 $\times {10}^{-5}$	3.9901 $\times {10}^{-7}$	1.2899 $\times {10}^{-7}$	8.9212 $\times {10}^{-10}$
4	1.3316 $\times {10}^{-5}$	4.0220 $\times {10}^{-8}$	4.1812 $\times {10}^{-8}$	3.3067 $\times {10}^{-9}$
14/3	4.6334 $\times {10}^{-5}$	5.0262 $\times {10}^{-7}$	7.2531 $\times {10}^{-8}$	5.6205 $\times {10}^{-9}$
16/3	7.0563 $\times {10}^{-6}$	6.0077 $\times {10}^{-7}$	8.5102 $\times {10}^{-8}$	2.1541 $\times {10}^{-9}$
6	2.3220 $\times {10}^{-5}$	5.8303 $\times {10}^{-7}$	5.7180 $\times {10}^{-8}$	2.8296 $\times {10}^{-10}$
20/3	5.3760 $\times {10}^{-7}$	1.0743 $\times {10}^{-6}$	3.0893 $\times {10}^{-8}$	2.9507 $\times {10}^{-9}$
22/3	1.9927 $\times {10}^{-5}$	9.4022 $\times {10}^{-7}$	1.3341 $\times {10}^{-8}$	1.7534 $\times {10}^{-9}$
8	6.2533 $\times {10}^{-6}$	7.1701 $\times {10}^{-7}$	4.5568 $\times {10}^{-8}$	1.3869 $\times {10}^{-9}$
26/3	1.4354 $\times {10}^{-5}$	5.8598 $\times {10}^{-7}$	4.7210 $\times {10}^{-8}$	1.9169 $\times {10}^{-9}$
28/3	3.9143 $\times {10}^{-7}$	5.9701 $\times {10}^{-7}$	5.4250 $\times {10}^{-9}$	2.9030 $\times {10}^{-9}$
10	6.2781 $\times {10}^{-6}$	3.4265 $\times {10}^{-7}$	2.4479 $\times {10}^{-8}$	1.1170 $\times {10}^{-9}$

and Figure 1 are examined the difference ${ϵ}_w=|w_{numerical}-{\overline{w}}_{OPIM}|$ for different values of the $N_{max}$ index.

From the Table 1 it results that with the increasing of $N_{max}$ index, the magnitude of the difference ${ϵ}_w=|w_{numerical}-{\overline{w}}_{OPIM}|$ decreases to ${10}^{-9}$. The optimal value for $N_{max}$ index towards to $N_{max}=28$.

The OPIM scheme is validated in the first case for the initial conditions: $x_0=0.25$, $y_0=0.75$, $z_0=1.05$. The profiles of the functions ${\overline{w}}_{OPIM}\in {\mathcal{SA}}_1$, $({\overline{x}}_{OPIM},{\overline{y}}_{OPIM},{\overline{z}}_{OPIM})\in {\mathcal{CFS}}_1$ are depicted in the Figure 2 and Figure 3.

The comparison between the functions: ${\overline{w}}_{OPIM}\in {\mathcal{SA}}_1$, $({\overline{x}}_{OPIM},{\overline{y}}_{OPIM},{\overline{z}}_{OPIM})\in {\mathcal{CFS}}_1$ and corresponding numerical ones, respectively, is shown in

Table 2. The numerical values of the function ${\overline{w}}_{OPIM}\in {\mathcal{SA}}_1$ of the Equation (6) using Equations (27) and (28) and corresponding numerical ones for initial data: $x_0=0.45$, $y_0=0.25$, $z_0=0.75$, physical parameters: $k=0.85$, $l=0.5$ and index $N_{max}=23$; (absolute errors: ${ϵ}_w=|w_{numerical}-{\overline{w}}_{OPIM}|$).

t	${\mathit{w}}_{\mathit{numerical}}$	${\overline{\mathit{w}}}_{\mathit{OPIM}}$	${\mathit{ϵ}}_{\mathit{w}}$
0	0.2893504211	0.2893504211	1.8596 $\times {10}^{-14}$
2/3	0.1032910999	0.1032911813	8.1428 $\times {10}^{-8}$
4/3	0.1011201382	0.1011200645	7.3662 $\times {10}^{-8}$
2	0.2831897278	0.2831896712	5.6622 $\times {10}^{-8}$
8/3	0.6273137782	0.6273138498	7.1642 $\times {10}^{-8}$
10/3	1.0843038063	1.0843039353	1.2899 $\times {10}^{-7}$
4	1.4174758199	1.4174758617	4.1812 $\times {10}^{-8}$
14/3	1.2690283219	1.2690282493	7.2531 $\times {10}^{-8}$
16/3	0.8179409154	0.8179408303	8.5102 $\times {10}^{-8}$
6	0.4137082021	0.4137081449	5.7180 $\times {10}^{-8}$
20/3	0.1584184708	0.1584184399	3.0893 $\times {10}^{-8}$
22/3	0.0789735113	0.0789735246	1.3341 $\times {10}^{-8}$
8	0.1875346655	0.1875347111	4.5568 $\times {10}^{-8}$
26/3	0.4682902423	0.4682902896	4.7210 $\times {10}^{-8}$
28/3	0.8912895213	0.8912895158	5.4250 $\times {10}^{-9}$
10	1.3240366947	1.3240366703	2.4479 $\times {10}^{-8}$

, Figure 2 and Figure 3, respectively.

4.2. Case 2: $k<0$, $l>0$

Figure 4 shows the variation of the difference ${ϵ}_w=|w_{numerical}-{\overline{w}}_{OPIM}|$ in the second case for initial data $x_0=0.55$, $y_0=0.35$, $z_0=0.25$, physical parameters $k=-0.65$, $l=0.75$ and index $N_{max}=23$.

Also, in this case the accuracy of the applied method is highlighted by comparison between the functions: ${\overline{w}}_{OPIM}\in {\mathcal{SA}}_2$, respectively $({\overline{x}}_{OPIM},{\overline{y}}_{OPIM},{\overline{z}}_{OPIM})\in {\mathcal{CFS}}_2$ and corresponding numerical ones in

Table 3. The numerical values of the function ${\overline{w}}_{OPIM}\in {\mathcal{SA}}_2$ using Equations (38) and (39) and corresponding numerical ones for initial data: $x_0=0.55$, $y_0=0.35$, $z_0=0.25$, physical parameters: $k=-0.65$, $l=0.75$ and index $N_{max}=23$; (absolute errors: ${ϵ}_w=|w_{numerical}-{\overline{w}}_{OPIM}|$).

t	${\mathit{w}}_{\mathit{numerical}}$	${\overline{\mathit{w}}}_{\mathit{OPIM}}$	${\mathit{ϵ}}_{\mathit{w}}$
0	0.3401219466	0.3401219466	1.4332 $\times {10}^{-13}$
1	0.4736249701	0.4736249707	5.6057 $\times {10}^{-10}$
2	0.4663420042	0.4663420037	5.3576 $\times {10}^{-10}$
3	0.3295734711	0.3295734711	1.6075 $\times {10}^{-11}$
4	0.2607145843	0.2607145843	2.0405 $\times {10}^{-11}$
5	0.3489310602	0.3489310601	7.3811 $\times {10}^{-11}$
6	0.4788049550	0.4788049552	2.4489 $\times {10}^{-10}$
7	0.4598070583	0.4598070584	1.1693 $\times {10}^{-10}$
8	0.3213952917	0.3213952916	7.2704 $\times {10}^{-11}$
9	0.2618232075	0.2618232075	1.9846 $\times {10}^{-11}$
10	0.3579406175	0.3579406174	9.9838 $\times {10}^{-11}$

and Figure 5, respectively Figure 6.

4.3. Case 3: $k\ne 0$, $l<0$

In Figure 7 is presented the variation of the difference ${ϵ}_w=|w_{numerical}-{\overline{w}}_{OPIM}|$ in the last case for initial data $x_0=0.55$, $y_0=0.35$, $z_0=0.25$, physical parameters $k=-0.450$, $l=-0.25$ and index $N_{max}=23$.

A comparison between the functions: ${\overline{w}}_{OPIM}\in {\mathcal{SA}}_3$, $({\overline{x}}_{OPIM},{\overline{y}}_{OPIM},{\overline{z}}_{OPIM})\in {\mathcal{SA}}_3$ and corresponding numerical ones, respectively, is emphasized in

Table 4. The numerical values of the function ${\overline{w}}_{OPIM}\in {\mathcal{SA}}_3$ using Equations (47) and (48) and corresponding numerical ones for initial data: $x_0=0.55$, $y_0=0.35$, $z_0=0.25$, physical parameters: $k=-0.450$, $l=-0.25$ and index $N_{max}=23$; (absolute errors: ${ϵ}_w=|w_{numerical}-{\overline{w}}_{OPIM}|$).

t	${\mathit{w}}_{\mathit{numerical}}$	${\overline{\mathit{w}}}_{\mathit{OPIM}}$	${\mathit{ϵ}}_{\mathit{w}}$
0	1.2120256565	1.2120256564	2.8498 $\times {10}^{-11}$
1	1.6123394336	1.6123394313	2.3431 $\times {10}^{-9}$
2	2.1997295278	2.1997295281	2.3885 $\times {10}^{-10}$
3	2.7369569941	2.7369569938	3.4806 $\times {10}^{-10}$
4	3.0163903678	3.0163903678	2.8598 $\times {10}^{-11}$
5	2.9493294559	2.9493294579	2.0125 $\times {10}^{-9}$
6	2.5553295952	2.5553295977	2.5374 $\times {10}^{-9}$
7	1.9697421940	1.9697421948	7.7861 $\times {10}^{-10}$
8	1.4287007553	1.4287007563	9.8604 $\times {10}^{-10}$
9	1.1420051499	1.1420051443	5.5923 $\times {10}^{-9}$
10	1.2010666437	1.2010666513	7.6534 $\times {10}^{-9}$

and Figure 8, respectively Figure 9.

The variation of the difference ${ϵ}_w=|w_{numerical}-{\overline{w}}_{OPIM}|$ for initial data $x_0=0.85$, $y_0=0.95$, $z_0=0.35$, physical parameters $k=0.750$, $l=-0.45$ and index $N_{max}=16$ is graphically plotted in Figure 10.

From the last couple of table and figures (

Table 5. The function ${\overline{w}}_{OPIM}\in {\mathcal{SA}}_3$ using Equations (47) and (48) and corresponding numerical ones for initial data: $x_0=0.85$, $y_0=0.95$, $z_0=0.35$, physical parameters: $k=0.750$, $l=-0.45$ and index $N_{max}=16$; (absolute errors: ${ϵ}_w=|w_{numerical}-{\overline{w}}_{OPIM}|$).

t	${\mathit{w}}_{\mathit{numerical}}$	${\overline{\mathit{w}}}_{\mathit{OPIM}}$	${\mathit{ϵ}}_{\mathit{w}}$
0	0.4636476090	0.4636476077	1.2868 $\times {10}^{-9}$
9/5	0.8460278121	0.8460287376	9.2546 $\times {10}^{-7}$
18/5	0.9723999843	0.9724014226	1.4383 $\times {10}^{-6}$
27/5	1.0173122060	1.0173114776	7.2840 $\times {10}^{-7}$
36/5	1.0438455818	1.0438430202	2.5615 $\times {10}^{-6}$
9	1.0893262516	1.0893231020	3.1496 $\times {10}^{-6}$
54/5	1.2176298207	1.2176255150	4.3056 $\times {10}^{-6}$
63/5	1.6058674831	1.6058603571	7.1260 $\times {10}^{-6}$
72/5	2.7109554512	2.7109452066	1.0244 $\times {10}^{-5}$
81/5	4.8014656171	4.8014740981	8.4809 $\times {10}^{-6}$
18	6.3883807767	6.3883687842	1.1992 $\times {10}^{-5}$

and Figure 11, respectively Figure 12) is found a good agreement between the OPIM solutions: ${\overline{w}}_{OPIM}$, respectively ${\overline{x}}_{OPIM}$, ${\overline{y}}_{OPIM}$, ${\overline{z}}_{OPIM}$ and corresponding numerical ones.

The numerical values of the convergence-control parameters for ${\overline{w}}_{OPIM}$ for four values of the $N_{max}$ index are presented in details in Appendix A.

4.4. OPIM Solutions versus Iterative Solutions

This subsection is dedicated for a comparison between the obtained results by OPIM scheme (using one iteration) and corresponding results with the iterative method described in [30] (using seven iterations). In the

Table 6. Values of the approximate analytical solution ${\overline{x}}_{OPIM}\left(t\right)$ (A3), the iterative solution $x_{iter}\left(t\right)$ (51) and the corresponding numerical ones.

t	${\mathit{x}}_{\mathit{numerical}}$	${\overline{\mathit{x}}}_{\mathit{OPIM}}$	${\mathit{x}}_{\mathit{iter}}$
0	0.45	0.45	0.45
0.35	0.5605299381	0.5605299043	0.5605299603
0.7	0.6019399952	0.6019399315	0.6019389881
1.05	0.6097536372	0.6097534668	0.6097448748
1.4	0.5961041764	0.5961041112	0.5960756769
1.75	0.5416320669	0.5416320764	0.5416749615
2.1	0.4088241222	0.4088245497	0.4101268423
2.45	0.1751981381	0.1751989192	0.1832389471
2.8	−0.1336579636	−0.1336581146	−0.1101295639
3.15	−0.4483418224	−0.4483418574	−0.4177692516
3.5	−0.6986536891	−0.6986538977	−0.6351428426

and Figure 13, Figure 14 and Figure 15, respectively, can be seen the precision and efficiency of the OPIM method.

The steps of the iterative method [30] are presented below.

If the system (1) is integrated over the interval $[0,t]$, it results:

$$\begin{array}{c}{\displaystyle x\left(t\right)=x\left(0\right)+\underset{0}{\overset{t}{\int }}z\left(s\right)\left(k-y\left(s\right)\right)ds}\hfill \\ {\displaystyle y\left(t\right)=y\left(0\right)+\underset{0}{\overset{t}{\int }}z\left(s\right)\left(l+x\left(s\right)\right)ds}\hfill \\ {\displaystyle z\left(t\right)=z\left(0\right)+\underset{0}{\overset{t}{\int }}\left(-kx\left(s\right)-ly\left(s\right)\right)ds}\hfill \end{array}.$$

(49)

The iterative algorithm leads to:

$$\begin{array}{c}{\displaystyle x_0\left(t\right)=x\left(0\right),x_1\left(t\right)=N_1(x_0,y_0,z_0)=\underset{0}{\overset{t}{\int }}z\left(0\right)\left(k-y\left(0\right)\right)ds,}\hfill \\ {\displaystyle y_0\left(t\right)=y\left(0\right),y_1\left(t\right)=N_2(x_0,y_0,z_0)=\underset{0}{\overset{t}{\int }}z\left(0\right)\left(l+x\left(0\right)\right)ds,}\hfill \\ {\displaystyle z_0\left(t\right)=z\left(0\right),z_1\left(t\right)=N_3(x_0,y_0,z_0)=\underset{0}{\overset{t}{\int }}\left(-kx\left(0\right)-ly\left(0\right)\right)ds,}\hfill \\ \cdots \hfill \\ x_m\left(t\right)=N_1\left({\displaystyle \sum _{i=0}^{m-1}{x}_i,\sum _{i=0}^{m-1}{y}_i,\sum _{i=0}^{m-1}{z}_i}\right)-N_1\left({\displaystyle \sum _{i=0}^{m-2}{x}_i,\sum _{i=0}^{m-2}{y}_i,\sum _{i=0}^{m-2}{z}_i}\right),\hfill \\ y_m\left(t\right)=N_2\left({\displaystyle \sum _{i=0}^{m-1}{x}_i,\sum _{i=0}^{m-1}{y}_i,\sum _{i=0}^{m-1}{z}_i}\right)-N_2\left({\displaystyle \sum _{i=0}^{m-2}{x}_i,\sum _{i=0}^{m-2}{y}_i,\sum _{i=0}^{m-2}{z}_i}\right),\hfill \\ z_m\left(t\right)=N_3\left({\displaystyle \sum _{i=0}^{m-1}{x}_i,\sum _{i=0}^{m-1}{y}_i,\sum _{i=0}^{m-1}{z}_i}\right)-N_3\left({\displaystyle \sum _{i=0}^{m-2}{x}_i,\sum _{i=0}^{m-2}{y}_i,\sum _{i=0}^{m-2}{z}_i}\right),\hfill \\ m\ge 2.\hfill \end{array}$$

(50)

The solutions of the Equations (1), using the iterative algorithm, can be written as:

$${\displaystyle x_{iter}\left(t\right)=\sum _{m=0}^{\infty }{x}_m\left(t\right),y_{iter}\left(t\right)=\sum _{m=0}^{\infty }{y}_m\left(t\right),z_{iter}\left(t\right)=\sum _{m=0}^{\infty }{z}_m\left(t\right),}$$

The iterative solutions $x_{iter}\left(t\right)$, after seven iterations and considering the initial conditions: $x\left(0\right)=0.45$, $y\left(0\right)=0.25$, $z\left(0\right)=0.75$ (presented in the Table 5), taking into account the algorithm (50), become:

$${\displaystyle \begin{array}{c}{\displaystyle x_{iter}\left(t\right)=\sum _{m=0}^7{x}_m\left(t\right)=0.45+0.4499999999t-0.4194375t^2+0.0647343750t^3+}\hfill \\ +0.1321203906t^4-0.0596508925t^5-0.0252664417t^6+0.0199209990t^7+\hfill \\ +0.0009944516t^8-0.0050189981t^9+0.0010096911t^{10}+0.0008782764t^{11}-\hfill \\ -0.0004015577t^{12}-0.0000830888t^{13}+0.0000926545t^{14}-5.042657·{10}^{-6}{t}^{15}-\hfill \\ -0.0000147286t^{16}+3.825744·{10}^{-6}{t}^{17}+1.521174·{10}^{-6}{t}^{18}-8.658214·{10}^{-7}{t}^{19}-\hfill \\ -4.960023·{10}^{-8}{t}^{20}+1.266344·{10}^{-7}{t}^{21}-1.683638·{10}^{-8}{t}^{22}-1.258976·{10}^{-8}{t}^{23}+\hfill \\ +\mathcal{O}\left({10}^{-9}\right),\hfill \\ y_{iter}\left(t\right)={\displaystyle \sum _{m=0}^7}{y}_m\left(t\right)=0.25+0.7124999999t-0.0723125000t^2-0.2979531249t^3+\hfill \\ +0.0548861588t^4+0.0604765449t^5-0.0344306378t^6-0.0068848995t^7+\hfill \\ +0.0108222194t^8-0.0007986684t^9-0.0022466787t^{10}+0.0007046456t^{11}+\hfill \\ +0.0003046745t^{12}-0.0002106980t^{13}-0.0000119912t^{14}+0.0000406429t^{15}-\hfill \\ -6.557029·{10}^{-6}{t}^{16}-5.337096·{10}^{-6}{t}^{17}+2.135663·{10}^{-6}{t}^{18}+3.785304·{10}^{-7}{t}^{19}-\hfill \\ -3.834523·{10}^{-7}{t}^{20}+2.179438·{10}^{-8}{t}^{21}+4.648639·{10}^{-8}{t}^{22}-1.140559·{10}^{-8}{t}^{23}+\hfill \\ +\mathcal{O}\left({10}^{-9}\right),\hfill \\ z_{iter}\left(t\right)={\displaystyle \sum _{m=0}^7}{z}_m\left(t\right)=0.75-0.5075000000t-0.3693749999t^2+0.1308927083t^3+\hfill \\ +0.0234880859t^4-0.0279490822t^5+0.0034108310t^6+0.0055273991t^7-\hfill \\ -0.0017171265t^8-0.0006315664t^9+0.0004262928t^{10}+7.531273·{10}^{-6}{t}^{11}-\hfill \\ -0.0000680699t^{12}+0.0000124592t^{13}+6.416595·{10}^{-6}{t}^{14}-2.766568·{10}^{-6}{t}^{15}-\hfill \\ -1.797293·{10}^{-7}{t}^{16}+3.239642·{10}^{-7}{t}^{17}-4.190182·{10}^{-8}{t}^{18}-2.083248·{10}^{-8}{t}^{19}+\hfill \\ +\mathcal{O}\left({10}^{-9}\right).\hfill \end{array}}$$

(51)

A comparison between the OPIM solutions ${\overline{x}}_{OPIM}$, ${\overline{y}}_{OPIM}$, ${\overline{z}}_{OPIM}$, the corresponding iterative solutions given in Equation (51) and the corresponding numerical ones are represented in Figure 13, Figure 14 and Figure 15 and in Table 5. The accuracy of the OPIM approach is highlighted by this analysis. This comparative analysis highlights the efficiency and the accuracy of the OPIM method using only one iteration.

5. Conclusions

This work proposed an analytical approach: OPIM to solve the second-order nonlinear differential equations using just one iteration. It was chosen as a nonlinear dynamical system to obtain the analytical approximate solutions, the ball–plate problem, that possess a Hamilton–Poisson structure. By comparison, the analytical results with the corresponding numerical ones yields a good agreement between them. The comparison proves the accuracy of the applied method, the obtained solutions are approaching the exact solution. The advantages of the OPIM method comparative with other analytical methods is synthesized by: the writing of the solutions in an effective form, the efficiency, the convergence control (in the sense that the residual functions is smaller than 1) and the non-existence of small parameters.

The achievement of these results encouraged the study of dynamical systems with similar properties. The selection of various control algorithms could be one of the directions for future study.