Sage Journals: Discover world-class research

Abstract

A nonlinear vibration system arises in physics. Besides its mathematical model, it is of great importance to have an accurate and reliable solution to the system. Though there are many analytical methods, such as the variational iteration method and the homotopy perturbation method, numerical approaches are rare. This paper suggests the barycentric interpolation collocation method to solve nonlinear oscillators. The Duffing equation is adopted as an example to elucidate the solution process. Some numerical examples are studied to demonstrate the accuracy of the present method. Results obtained by the method indicate the method is simple and effective.

Keywords

Duffing equation vibration system direct linearized iterative barycentric interpolation collocation method numerical solution

Introduction

The nonlinear vibration system has wide applications in physics. There are some valuable efforts that focus on finding analytical methods for solving the nonlinear vibration system. These analytical methods are the variational iteration method (VIM), the homotopy perturbation method (HPM), the Adomian decomposition method and so on. The VIM and the HPM have proved to be powerful mathematical tools for solving various kinds of nonlinear problems.^1–8 The barycentric interpolation collocation method (BICM)^9–11 is a high precision numerical method. Some authors^12,13 have used the BICM to solve some linear and nonlinear high-dimensional Fredholm integral equations,¹² nonlinear parabolic partial differential equations.^13,14 In this paper, the BICM is used to solve nonlinear oscillators and coupled Duffing system is adopted as an example to elucidate the solution process.

We consider the following coupled nonlinear Duffing system¹⁵

\begin{array}{l} {\begin{array}{l} \frac{d^{2} x}{d t^{2}} = - ε (K + 1) \frac{d x}{d t} - x - x^{3} + ε K \frac{d x}{d t} + B \cos ω t \\ \frac{d^{2} y}{d t^{2}} = - ε (K + 1) \frac{d y}{d t} - y - y^{3} + ε K \frac{d y}{d t} + B \cos ω t \end{array} \end{array}

(1)

with the initial conditions

x (0) = c_{1}, y (0) = c_{2}, x' (0) = c_{3}, y' (0) = c_{4}

(2)

where

K, ε, B, ω, c_{1}, c_{2}, c_{3}, c_{4}

are real constants,

0 < t < T

x (t), y (t)

are unknown functions.

The numerical solution of (1)

The interval $[0, T]$ is divided into M different nodes: $0 = t_{1} < t_{2} < \dots < t_{M} = T$ . For $i = 1, 2, \dots, N$ , we give two initial hypothesis functions $x_{0} (t), y_{0} (t)$ , and construct following a linear iterative format of the model (1)

\begin{array}{l} {\begin{array}{l} \frac{d^{2} x_{n}}{d t^{2}} + ε (K + 1) \frac{d x_{n}}{d t} + x_{n} - ε K \frac{d y_{n}}{d t} = - x_{n - 1}^{3} + B \cos ω t \\ \frac{d^{2} y_{n}}{d t^{2}} + ε (K + 1) \frac{d y_{n}}{d t} + y_{n} - ε K \frac{d x_{n}}{d t} = - y_{n - 1}^{3} + B \cos ω t \end{array} \end{array} n = 1, 2, 3, \dots

(3)

Let

x_{n} (t) = \sum_{j = 1}^{M} ξ_{j} (t) x_{n} (t_{j}), y_{n} (t) = \sum_{j = 1}^{M} ξ_{j} (t) y_{n} (t_{j})

(4)

where

ξ_{j} (t)

is the barycentric Lagrange interpolation primary function, ω_j is the center of gravity interpolation weight.

ξ_{j} (t) = \frac{\frac{ω_{j}}{t - t_{j}}}{\sum_{k = 1}^{M} \frac{ω_{k}}{t - t_{k}}}, ω_{j} = \frac{1}{\prod_{i = 1, j \neq k}^{M} (t_{i} - t_{j})}

(5)

So, when $t = t_{i}, i = 1, 2, \dots M$ , the format (3) can be written in the following matrix form

[\begin{matrix} D ″ + ε (K + 1) D' + I & - ε K D' \\ - ε K D' & D ″ + ε (K + 1) D' + I \end{matrix}] [\begin{matrix} x_{n} \\ y_{n} \end{matrix}] = [\begin{matrix} - x_{n - 1}^{3} + B \cos ω t \\ - y_{n - 1}^{3} + B \cos ω t \end{matrix}]

(6)

where the matrix

D' = {(\frac{d ξ_{j} (t)}{d t} |_{t = t_{i}})}_{i, j = 1, 2, \dots M}

, and

D ″ = {(\frac{d^{2} ξ_{j} (t)}{d t^{2}} |_{t = t_{i}})}_{i, j = 1, 2, \dots M,}

are the M order matrices. I is the M order unit matrix. The vector

\begin{array}{l} [x_{n}, y_{n}] = [x_{n} (t_{1}), x_{n} (t_{2}), \dots, x_{n} (t_{M}), y_{n} (t_{1}), y_{n} (t_{2}), \dots, y_{n} (t_{M})], \\ [- x_{n - 1}^{3} + B \cos ω t, - y_{n - 1}^{3} + B \cos ω t] = [- x_{n - 1}^{3} (t_{1}) + B \cos ω t_{1}, - x_{n - 1}^{3} (t_{2}) \\ + B \cos ω t_{2}, \dots, - x_{n - 1}^{3} (t_{M}) + B \cos ω t_{M}, - y_{n - 1}^{3} (t_{1}) + B \cos ω t_{1}, - y_{n - 1}^{3} (t_{2}) \\ + B \cos ω t_{2}, \dots, - y_{n - 1}^{3} (t_{M}) + B \cos ω t_{M}] \end{array}

Next, we take the formula (4) into the initial conditions (2), we can get the following equation

\sum_{i = 1}^{M} ξ_{i} (0) x_{n} (t_{i}) = c_{1}, \sum_{i = 1}^{M} ξ_{i}^{'} (0) x_{n} (t_{i}) = c_{3}, \sum_{i = 1}^{M} ξ_{i} (0) y_{n} (t_{i}) = c_{2}, \sum_{i = 1}^{M} ξ_{i}^{'} (0) y_{n} (t_{i}) = c_{4}

(7)

The first 1 and the first 2M of equation (5) are replaced separately by equation (7) in turn. So, we can get $x_{n} (t), y_{n} (t)$ as approximate solution of (1) and (2).

Numerical experiments

In this section, some numerical examples are studied to demonstrate the accuracy of the present method. The examples are computed using MatlabR2017a.

Example 1 We consider the following the Duffing equation

\frac{d^{2} x}{d t^{2}} + ε \frac{d x}{d t} + x + b x^{3} = B \cos ω t, t \in [0, T]

(8)

with the following initial conditions

x (0) = 0, x' (0) = 1

(9)

In Figures 1 to 6, we give some numerical results of Example 1. The parameters of Figures 1 to 6 are listed in Table 1.

Table 1.

Parameters used in the numerical simulations for Example 1, $b = 1, ε = 0.18,$ $ω = 0.8$ .

Figure	t	B
Figure 1	[0,1]	18.5
Figure 2	[0,100]	18.5
Figure 3	[0,1]	23.5
Figure 4	[0,100]	23.5
Figure 5	[0,1]	26.7
Figure 6	[0,100]	26.7

Figure 1.

Numerical results for Example 1 at B = 18.5, $0 \leq t \leq 1$ .

Figure 2.

Numerical results for Example 1 at B = 18.5, $0 \leq t \leq 100$ .

Figure 3.

Numerical results for Example 1 at B = 23.5, $0 \leq t \leq 1$ .

Figure 4.

Numerical results for Example 1 at B = 23.5, $0 \leq t \leq 100$ .

Figure 5.

Numerical results for Example 1 at B = 26.7, $0 \leq t \leq 1$ .

Figure 6.

Numerical results for Example 1 at B = 26.7, $0 \leq t \leq 100$ .

Example 2 We consider the model (1) with the initial conditions $x (0) = 0, x' (0) = 0, y (0) = 1, y' (0) = 0.$

In Figures 7 to 12, we give some numerical results of Example 2. The parameters of Figures 7 to 12 are listed in Table 2.

Table 2.

Parameters used in the numerical simulations for Example 2, $0 \leq t \leq 100$ .

Figure	ε	ω	B	K
Figure 7	0.18	0.8	23.5	5
Figure 8	0.18	0.8	23.5	200
Figure 9	5.18	0.8	23.5	5
Figure 10	0.18	3.4	23.5	5
Figure 11	0.18	0.8	5	5
Figure 12	0.18	3.4	5	5

Figure 7.

Numerical results for Example 2 at $B = 23.5, K = 5, 0 \leq t \leq 100$ .

Figure 8.

Numerical results for Example 2 at $B = 23.5, K = 200, 0 \leq t \leq 100$ .

Figure 9.

Numerical results for Example 2 at $B = 23.5, K = 5, ε = 5.18, 0 \leq t \leq 100$ .

Figure 10.

Numerical results for Example 2 at $B = 23.5, K = 5, w = 3.4, 0 \leq t \leq 100$ .

Figure 11.

Numerical results for Example 2 at $B = 5, K = 5, 0 \leq t \leq 100$ .

Figure 12.

Numerical results for Example 2 at $B = 5, K = 5, w = 3.4, 0 \leq t \leq 100$ .

From the above Examples, we can see that the BICM is an effective method. Next, we use the BICM to solve the following nonlinear epidemic model.

Example 3 We consider the following nonlinear epidemic model

\begin{array}{l} {\begin{array}{l} \frac{d x}{d t} = - 0.1 x (t) y (t) \\ \frac{d y}{d t} = 0.1 x (t) y (t) - 0.2 y (t), 0 < t < T \\ \frac{d z}{d t} = 0.2 y (t) \end{array} \end{array}

(10)

with the following initial conditions

x (0) = c_{1}, y (0) = c_{2}, z (0) = c_{3}

(11)

x(t) susceptible population: those so far uninfected and therefore liable to infection; y(t) infective population: those who have the disease and are still at large; z(t) isolated population, or who have recovered and are therefore immune. The numerical results of Example 3 are shown in Figure 13(a) to (c). The values of c_1, c₂ and c₃ are shown in Table 3. From Figure 13, in a short time t, we can know that:

Figure 13.

Changes in population composition values affect the epidemic model for Example 3 (a), Example 3 (b), and Example 3 (c).

Table 3.

Parameter selection for Example 3.

c ₁	c ₂	c ₃
15	15	15
0	15	10
15	20	10

When we take the same value for susceptible population, infective population and isolated population, the rate of decrease of susceptible population is fast, and infective population is always more than isolated population.

When we take the value of susceptible population is zero, the infective population is the largest, the infective population is inversely proportional to the isolated population.

When we take the value of infective population is the largest, the last infective population and the susceptible population are stable.

Conclusions and remarks

In this paper, the Duffing systems have been solved by the BICM. From some figures, we can see that the present method is feasible in solving this model. Meanwhile, these numerical experiments illustrate that the numerical results of the present method are the same as the experimental results of Tafo Wembe and Yamapi¹⁵ Compared with Tafo Wembe and Yamapi,¹⁵ we get better numerical results.

All computations are performed by the MatlabR2014a and MatlabR2017a software package.

Data availability

The data used to support the findings of this study are available from the corresponding author upon request.

Footnotes

Acknowledgements

The authors thank the reviewers for their valuable suggestions, which greatly improved the quality of the paper.

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This paper is supported by the Natural Science Foundation of Inner Mongolia (2017MS0103), Inner Mongolia Maker Collaborative Innovation Center of Jining Normal University and the Numerical Analysis of Graduate Course Construction Project of Inner Mongolia University of Technology (Grant Number KC2014001) and the National Natural Science Foundation of China (11361037).

References

JH.

Variational iteration method for delay differential equations. Commun Nonlinear Sci Numer Simul 1997; 2: 235–236.

JH.

Approximate solution of nonlinear differential equations with convolution product nonlinearities. Comput Methods Appl Mech Eng 1998; 167: 69–73.

JH.

Homotopy perturbation method: a new nonlinear analytical technique. Appl Math Comput 2003; 135: 73–79.

JH.

A coupling method of homotopy technique and perturbation technique for nonlinear problems. Int J Non-Linear Mech 2000; 35: 37–43.

JH.

Variational iteration method a kind of non-linear analytical technique: some examples. Int J Non-Linear Mech 1999; 34: 699–708.

JH.

Variational iteration method for autonomous ordinary differential systems. Appl Math Comput 2000; 114: 115–123.

XH.

Construction of solitary solution and compacton-like solution by variational iteration method. Chaos Solitons Fract 2006; 29: 108–113.

JH.

Some asymptotic methods for strongly nonlinear equations. Int J Mod Phys B 2006; 20: 1141–1199.

Wang

ZQ.

Barycentric interpolation collocation method for nonlinear problems. Beijing, China: National Defense Industry Press, 2015.

10.

Wang

ZQ.

High-precision non-grid center of gravity interpolation collocation method: algorithm, program and engineering application.Beijng, China: Science Press, 2012.

11.

Wang

Tian

ZY.

Numerical method for singularly perturbed delay parabolic partial differential equations. Therm Sci 2017; 21: 1595–1599.

12.

Liu

Huang

Pan

, et al. Barycentric interpolation collocation methods for solving linear and nonlinear high-dimensional Fredholm integral equations. Comput Math Appl 2018; 327: 141–154.

13.

Liu

Wang

SG.

Barycentric interpolation collocation method for solving the coupled viscous Burgers’ equations. Int J Comput Math 2018; 95: 2162–2173.

14.

Luo

Huang

, et al. Barycentric rational collocation methods for a class of nonlinear parabolic partial differential equations. Appl Math Lett 2017; 68: 13–19.

15.

Tafo Wembe

Yamapi

Chaos synchronization of resistively coupled Duffing systems: numerical and experimental investigations. Commun Nonlinear Sci Numer Simul 2009; 14: 1439–1453.

The barycentric interpolation collocation method for a class of nonlinear vibration systems

Abstract

Keywords

Introduction

The numerical solution of (1)

Numerical experiments

Conclusions and remarks

Data availability

Footnotes

Acknowledgements

Declaration of conflicting interests

Funding

References