Perturbation Methods to Analysis of Thermal, Fluid Flow and Dynamics Behaviors of Engineering Systems

Abstract

This chapter presents the applications of perturbation methods such as regular and homotopy perturbation methods to thermal, fluid flow and dynamic behaviors of engineering systems. The first example shows the utilization of regular perturbation method to thermal analysis of convective-radiative fin with end cooling and thermal contact resistance. The second example is concerned with the application of homotopy perturbation method to squeezing flow and heat transfer of Casson nanofluid between two parallel plates embedded in a porous medium under the influences of slip, Lorentz force, viscous dissipation and thermal radiation. Additionally, the dynamic behavior of piezoelectric nanobeam embedded in linear and nonlinear elastic foundations operating in a thermal-magnetic environment is analyzed using homotopy perturbation method which is presented in the third example. It is believed that the presentation in this chapter will enhance the understanding of these methods for the real world applications.

1. Introduction

The descriptions of the behaviors of the real world phenomena and systems through the use of mathematical models often involve developments of nonlinear equations which are difficult to solve exactly and analytically. Consequently, recourse is always made to numerical methods as alternative methods in solving the nonlinear equations. However, the developments of analytical solutions are obviously still very important. Analytical solutions for specified problems are also essential and required to show the direct relationship between the models parameters. When analytical solutions are available, they provide good insights into the significance of various system parameters affecting the phenomena. Such solutions provide continuous physical insights than pure numerical or computation methods. Indisputably, analytical solutions are convenient for parametric studies, accounting for the physics of the problem and appear more appealing than the numerical solutions. Also, they help in reducing the computation and simulation costs as well as the task involved in the analysis of real-life problems.

Although, there is no general exact analytical method to solve all nonlinear problems, over the years, the nonlinear problems have been solved using different approximate analytical methods such as regular perturbation, singular perturbation method, homotopy perturbation method, homotopy analysis method, methods of weighted residual, variational iterative method, differential transformation method, variation parameter method, Adomian decomposition method, etc. The non-perturbative approximate analytic methods present explicit approximate analytical solutions which often involve complex mathematical analysis leading to analytic expressions involving large number terms. Furthermore, the methods are inherently with high computational cost and time accompanied with the requirement of high skills in mathematics. Moreover, in practice, analytical solutions with large number of terms and conditional statements for the solutions are not convenient for use by designers and engineers. Also, in these methods, there are always search for particular value(s) that will satisfy the end boundary condition(s). This always necessitates the use of software and such could result in additional computational cost in the generation of solution to the problem. Also, the quests involve applications of numerical schemes to determine the required value(s) that will satisfy the end boundary condition(s). This fact renders most of the approximate analytical methods to be taken as more of semi-analytical methods than total approximate analytical methods. Moreover, these methods have their own operational restrictions that severely narrow their functioning domain and when they are routinely implemented, they can sometimes lead to erroneous results. Specifically, the transformation of the nonlinear equations and the development of equivalent recurrence equations for the nonlinear equations using differential transformation method proved somehow difficult in some nonlinear system such as in rational Duffing oscillator, irrational nonlinear Duffing oscillator, finite extensibility nonlinear oscillator. There is difficulty in the determination of Adomian polynomials for the application of Adomian decomposition method for nonlinear problems. There are lack of rigorous theories or proper guidance for choosing initial approximation, auxiliary linear operators, auxiliary functions, and auxiliary parameters in the use of homotopy analysis method. Therefore, the need for comparatively simple, flexible, generic and high accurate total approximate analytical solutions is well established. One of the techniques that can be applied for such quest is the perturbation method. Perturbation method, although comparably old, as a pioneer method for finding approximate analytical solutions to nonlinear problems, it offers an alternative approach to solving certain types of nonlinear problems. In the limit of small parameter, perturbation method is widely used for solving many heat transfer, vibration, fluid mechanics and solid mechanics problems. It is capable of solving nonlinear, inhomogeneous and multidimensional problems with reasonable high level of accuracy. The most significant efforts and applications of the method were focused on celestial mechanics, fluid mechanics, and aerodynamics. Although, the solutions reported for other sophisticated methods to difference problems have good accuracy, they are more complicated for applications than perturbation method. Therefore, over the years, the relative simplicity and high accuracy especially in the limit of small parameter have made perturbation method an interesting tool among the most frequently used approximate analytical methods. Although, the perturbation method provides in general, better results for small perturbation parameters, besides having a handy mathematical formulation, it has been shown to have a good accuracy, even for relatively large values of the perturbation parameter [1, 2, 3, 4, 5].

2. Example 1: regular perturbation method to thermal analysis of convective-radiative fin with end cooling and thermal contact resistance

Consider a convective-radiative fin of temperature-dependent thermal conductivity k(T), length L and thickness δ, exposed on both faces to a convective environment at temperature $T_{\infty }$ and a heat transfer co-efficient h subjected to magnetic field shown in Figure 1. The dimension x pertains to the length coordinate which has its origin at the tip of the fin and has a positive orientation from the fin tip to the fin base. In order to analyze the problem, the following assumptions are made. The following assumptions were made in the development of the model

1. The heat flow in the fin and its temperatures remain constant with time.

2. The temperature of the medium surrounding the fin is uniform.

3. The temperature of the base of the fin is uniform.

4. The fin thickness is small compared with its width and length, so that temperature gradients across the fin thickness and heat transfer from the edges of the fin is negligible compared with the heat leaving its lateral surface.

Figure 1.

(a) Schematic of the convective-radiative longitudinal straight fin with magnetic field. (b) Schematic of the longitudinal straight fin geometry showing thermal contact resistance and boundary conditions.

Applying thermal energy balance on the fin and using the above model assumptions, the following nonlinear thermal model is developed

$$\frac{d}{\mathit{dx}}\left[\left[1+\lambda \left(T-T_a\right)\right]\frac{\mathit{dT}}{\mathit{dx}}\right]-\frac{h}{k_a\delta }\left(T-T_a\right)-\frac{\mathit{\sigma \epsilon }}{k_a\delta }\left(T^4-T_a^4\right)-\frac{\sigma B_o^2{u}^2}{k_a{A}_{\mathit{cr}}}\left(T-T_a\right)=0$$E1

The boundary conditions are

$$x=0,-k\left(T\right)\frac{\partial T}{\partial x}=h_e\left(T-T_a\right)+\sigma \in \left(T^4-T_a^4\right)$$E2

$$x=L,-k\left(T\right)\frac{\partial T}{\partial x}=h_c\left(T_b-T\right)+\sigma \in \left(T^4-T_a^4\right)$$E3

Considering a case when a small temperature difference exists within the material during the heat flow. This actually necessitated the use of temperature-invariant physical and thermal properties of the fin. Also, it has been established that under such scenario, the term T⁴ can be expressed as a linear function of temperature. Therefore, we have

$$T^4=T_a^4+4T_a^3\left(T-T_a\right)+6T_a^2{\left(T-T_a\right)}^2+\dots \cong 4T_a^{3}T-3T_a^4$$E4

On substituting Eq. (4) into Eq. (1), one arrives arrived at

$$\frac{d}{\mathit{dx}}\left[\left[1+\lambda \left(T-T_{\infty }\right)\right]\frac{\mathit{dT}}{\mathit{dx}}\right]-\frac{h}{k_a\delta }\left(T-T_a\right)-\frac{4\mathit{\sigma \epsilon }{T}_a^3}{k_a\delta }\left(T-T_a\right)-\frac{\sigma B_o^2{u}^2}{k_a{A}_{\mathit{cr}}}\left(T-T_a\right)=0$$E5

The boundary conditions

$$x=0,-k\left(T\right)\frac{\partial T}{\partial x}=h_e\left(T-T_a\right)+4\sigma \in T_a^3\left(T-T_a\right)$$E6

$$x=L,-k\left(T\right)\frac{\partial T}{\partial x}=h_c\left(T_b-T\right)+4\sigma \in T_a^3\left(T-T_a\right)$$E7

On introducing the following dimensionless parameters in Eq. (8) into Eq. (5),

$$X=\frac{x}{L},\theta =\frac{T-T_a}{T_b-T_a},\mathit{Ra}=\frac{\mathit{gk\beta }\left(T_b-T_a\right)b}{\mathit{\alpha \nu }{k}_r},N=\frac{4{\sigma }_{\mathit{st}}{\mathit{bT}}_a^3}{k_a},\mathit{Ha}=\frac{\sigma B_0^2{u}^2}{k_a{A}_{\mathit{cr}}}.$$E8

$$\begin{array}{l}{\mathit{Bi}}_e=\frac{h_{e}b}{k_a},{\mathit{Bi}}_c=\frac{h_{c}b}{k_a},M^2=\frac{{\mathit{hb}}^2}{k_a\delta },\epsilon =\lambda \left(T_b-T_a\right){\mathit{Bi}}_{e,\mathit{eff}}=\frac{\left(h_e+\mathit{\sigma \epsilon }\right)b}{k_a},\\ {\mathit{Bi}}_{\mathit{\text{ceff}}}=\frac{\left(h_c+\mathit{\sigma \epsilon }\right)b}{k_a}\end{array}$$

The dimensionless form of the governing Eq. (5) is arrived at as

$$\frac{d}{\mathit{dX}}\left[\left(1+\mathit{\epsilon \theta }\right]\frac{\mathit{d\theta }}{\mathit{dX}}\right]-M^2\theta -\mathit{Nr\theta }-\mathit{Ha\theta }=0$$E9

On expanding Eq. (9), one has

$$\frac{d^2\theta }{{\mathit{dX}}^2}+\mathit{\epsilon \theta }\frac{d^2\theta }{{\mathit{dX}}^2}+\epsilon {\left(\frac{\mathit{d\theta }}{\mathit{dX}}\right)}^2-M^2\theta -\mathit{Nr\theta }-\mathit{Ha\theta }=0$$E10

The boundary conditions are

$$X=0,\left(1+\mathit{\epsilon \theta }\right)\frac{\mathit{d\theta }}{\mathit{dX}}=-{\mathit{Bi}}_{e,\mathit{eff}}\theta$$E11

$$X=1,\left(1+\mathit{\epsilon \theta }\right)\frac{\mathit{d\theta }}{\mathit{dX}}=-{\mathit{Bi}}_{c,\mathit{eff}}\left(1-\theta \right)$$E12

3. Method of solution using regular perturbation method

It is very difficult to develop closed-form solution for the above non-linear Eq. (10). Therefore, in this work, recourse is made to apply a relatively simple and accurate method approximate analytical method, the perturbation method. Perturbation theory is based on the fact that the equation(s) describing the phenomena or process under investigation contain(s) a small parameter (or several small parameters), explicitly or implicitly. Therefore, the perturbation method is applicable to very small magnitudes of ε where the nonlinearity is slightly effective. Although, it has been shown to have a good accuracy, even for relatively large values of the perturbation parameter, ε [1, 2].

In solving Eq. (10), one needs to expand the dimensionless temperature as

$$\theta ={\theta }_0+{\mathit{\epsilon \theta }}_1+{\epsilon }^2{\theta }_2+\dots$$E13

Substituting Eq. (13) into Eq. (10), up to first order approximate, we have

$$\begin{array}{l}\frac{d^2{\theta }_0}{{\mathit{dX}}^2}-\left(M^2+\mathit{Nr}+\mathit{Ha}\right){\theta }_0+\epsilon \left[\frac{d^2{\theta }_1}{{\mathit{dX}}^2}+{\theta }_0\frac{d^2{\theta }_0}{{\mathit{dX}}^2}+{\left(\frac{d{\theta }_0}{\mathit{dX}}\right)}^2-\left(M^2+\mathit{Nr}+\mathit{Ha}\right){\theta }_1\right]\\ +{\epsilon }^2\left[\frac{d^2{\theta }_2}{{\mathit{dX}}^2}+{\theta }_1\frac{d^2{\theta }_0}{{\mathit{dX}}^2}+{\theta }_0\frac{d^2{\theta }_1}{{\mathit{dX}}^2}+2\left(\frac{d{\theta }_1}{\mathit{dX}}\right)\left(\frac{d{\theta }_0}{\mathit{dX}}\right)-\left(M^2+\mathit{Nr}+\mathit{Ha}\right){\theta }_2\right]=0\end{array}$$E14

Leading order and first order equations with the appropriate boundary conditions are given as:

Leading order equation:

$$\frac{d^2{\theta }_o}{{\mathit{dx}}^2}-\left(M^2+\mathit{Nr}+\mathit{Ha}\right){\theta }_o=0$$E15

Subject to:

$$X=0,\frac{d{\theta }_o}{\mathit{dX}}=-{\mathit{Bi}}_{e,\mathit{eff}}{\theta }_o$$E16

$$X=1,\frac{d{\theta }_o}{\mathit{dX}}={\mathit{Bi}}_{c,\mathit{eff}}\left({\theta }_o-1\right)$$E17

First-order equation:

$$\frac{d^2{\theta }_1}{{\mathit{dX}}^2}-\left(M^2+\mathit{Nr}+\mathit{Ha}\right){\theta }_1=-{\left(\frac{d{\theta }_o}{\mathit{dX}}\right)}^2-{\theta }_0\frac{d^2{\theta }_o}{{\mathit{dX}}^2}$$E18

Subject to:

$$X=0,{\theta }_0\frac{d{\theta }_0}{\mathit{dX}}+\frac{d{\theta }_1}{\mathit{dX}}=-{\mathit{Bi}}_{e,\mathit{eff}}{\theta }_1$$E19

$$X=1,{\theta }_0\frac{d{\theta }_0}{\mathit{dX}}+\frac{d{\theta }_1}{\mathit{dX}}={\mathit{Bi}}_{\mathrm{c},\mathit{eff}}{\theta }_1$$E20

Second-order equation

$$\frac{d^2{\theta }_2}{{\mathit{dX}}^2}-\left(M^2+\mathit{Nr}+\mathit{Ha}\right){\theta }_2=-{\theta }_1\frac{d^2{\theta }_0}{{\mathit{dX}}^2}-{\theta }_0\frac{d^2{\theta }_1}{{\mathit{dX}}^2}-2\left(\frac{d{\theta }_1}{\mathit{dX}}\right)\left(\frac{d{\theta }_0}{\mathit{dX}}\right)$$E21

The boundary conditions

$$X=0,{\theta }_1\frac{d{\theta }_o}{\mathit{dX}}+{\theta }_0\frac{d{\theta }_1}{\mathit{dX}}+\frac{d{\theta }_2}{\mathit{dX}}=-{\mathit{Bi}}_{e,\mathit{eff}}{\theta }_2$$E22

$$X=1,{\theta }_1\frac{d{\theta }_0}{\mathit{dX}}+{\theta }_0\frac{d{\theta }_1}{\mathit{dX}}+\frac{d{\theta }_2}{\mathit{dX}}={\mathit{Bi}}_{\mathrm{c},\mathit{eff}}{\theta }_2$$E23

It can be shown from Eq. (15), (18) and (21) with the corresponding boundary conditions of Eqs. (16), (19) and (22) that the:

Leading order solution for θ_o is

$${\theta }_o=\frac{{\mathit{Bi}}_c\left\{\sqrt{\left(M^2+\mathit{Nr}+\mathit{Ha}\right)}\mathit{cosh}\left(\sqrt{\left(M^2+\mathit{Nr}+\mathit{Ha}\right)}X\right)-{\mathit{Bi}}_e\mathrm{s}\mathit{inh}\left(\sqrt{\left(M^2+\mathit{Nr}+\mathit{Ha}\right)}X\right)\right\}}{\left\{\begin{array}{l}{\mathit{Bi}}_c\left\{\left(\sqrt{\left(M^2+\mathit{Nr}+\mathit{Ha}\right)}\right)\mathit{cosh}\left(\sqrt{\left(M^2+\mathit{Nr}+\mathit{Ha}\right)}\right)-{\mathit{Bi}}_e\mathit{sinh}\left(\sqrt{\left(M^2+\mathit{Nr}+\mathit{Ha}\right)}\right)\right\}\\ +\sqrt{\left(M^2+\mathit{Nr}+\mathit{Ha}\right)}\left\{{\mathit{Bi}}_e\mathit{cosh}\left(\sqrt{\left(M^2+\mathit{Nr}+\mathit{Ha}\right)}\right)-\left(\sqrt{\left(M^2+\mathit{Nr}+\mathit{Ha}\right)}\right)\mathit{sinh}\left(\sqrt{\left(M^2+\mathit{Nr}+\mathit{Ha}\right)}\right)\right\}\end{array}\right\}}$$E24

While the first order solution θ₁ is

E25

The second-order solution θ₂ is too huge to be included in the manuscript.

On substituting Eqs. (24) and (25) into Eq. (13) up to the first order (i.e. neglecting the higher orders), one arrives at

E26

4. Example 2: homotopy perturbation method to analysis of squeezing flow and heat transfer of Casson nanofluid between two parallel plates embedded in a porous medium under the influences of slip, Lorentz force, viscous dissipation and thermal radiation

Consider a Casson nanofluid flowing between two parallel plates placed at time-variant distance and under the influence of magnetic field as shown in the Figure 2. It is assumed that the flow of the nanofluid is laminar, stable, incompressible, isothermal, non-reacting chemically, the nanoparticles and base fluid are in thermal equilibrium and the physical properties are constant. The fluid conducts electrical energy as it flows unsteadily under magnetic force field. The fluid structure is everywhere in thermodynamic equilibrium and the plate is maintained at constant temperature.

Figure 2.

Model diagram of MHD squeezing flow of nanofluid between two parallel plates embedded in a porous medium.

Following the assumptions, the governing equations for the flow are given as

$$\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}=0$$E27

$${\rho }_{\mathit{nf}}\left(\frac{\partial u}{\partial t}+u\frac{\partial u}{\partial x}+v\frac{\partial u}{\partial y}\right)=-\frac{\partial p}{\partial x}+{\mu }_{\mathit{nf}}\left(1+\frac{1}{\beta }\right)\left(2\frac{{\partial }^{2}u}{\partial x^2}+\frac{{\partial }^{2}u}{\partial x\partial y}+\frac{{\partial }^{2}u}{\partial y^2}\right)-\sigma B_o^{2}u-\frac{{\mu }_{\mathit{nf}}u}{K_p}$$E28

$${\rho }_{\mathit{nf}}\left(\frac{\partial v}{\partial t}+u\frac{\partial v}{\partial x}+v\frac{\partial v}{\partial y}\right)=-\frac{\partial p}{\partial y}+{\mu }_{\mathit{nf}}\left(1+\frac{1}{\beta }\right)\left(2\frac{\partial v}{\partial x^2}+\frac{{\partial }^{2}v}{\partial x\partial y}+\frac{{\partial }^{2}v}{\partial y^2}\right)-\frac{{\mu }_{\mathit{nf}}v}{K_p}$$E29

$$\begin{array}{l}\frac{\partial T}{\partial t}+u\frac{\partial T}{\partial x}+u\frac{\partial T}{\partial y}=\frac{k_{\mathit{nf}}}{{\left(\rho C_p\right)}_{\mathit{nf}}}\left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}\right)\\ +\frac{{\mu }_{\mathit{nf}}}{{\left(\rho C_p\right)}_{\mathit{nf}}}\left(1+\frac{1}{\beta }\right)\left(\begin{array}{l}2{\left(\frac{{\partial }^{2}u}{\partial x^2}\right)}^2+{\left(\frac{{\partial }^{2}u}{\partial y^2}+\frac{{\partial }^{2}v}{\partial x^2}\right)}^2\\ +2{\left(\frac{{\partial }^{2}v}{\partial y^2}\right)}^2\end{array}\right)-\frac{1}{{\left(\rho C_p\right)}_{\mathit{nf}}}\frac{\partial q_r}{\partial x}\end{array}$$E30

where

$${\rho }_{\mathit{nf}}={\rho }_f\left(1-\varphi \right)+{\rho }_s\varphi$$E31

$${\mu }_{\mathit{nf}}=\frac{{\mu }_f}{{\left(1-\varphi \right)}^{2.5}}$$E32

and the magnetic field parameter is given as

$$B\left(t\right)=\frac{B_0}{\sqrt{1-\mathit{\alpha t}}}$$E33

$${\sigma }_{\mathit{nf}}={\sigma }_f\left[1+\frac{3\left\{\frac{{\sigma }_s}{{\sigma }_f}-1\right\}\varphi }{\left\{\frac{{\sigma }_s}{{\sigma }_f}+2\right\}\varphi -\left\{\frac{{\sigma }_s}{{\sigma }_f}-1\right\}\varphi }\right],$$E34

$$k_{\mathit{nf}}=k_f\left[\frac{k_s+\left(m-1\right)k_f-\left(m-1\right)\left(k_f-k_s\right)\varphi }{k_s+\left(m-1\right)k_f+\left(k_f-k_s\right)\varphi }\right],$$E35

The Casson fluid parameter, $\beta ={\mu }_B\sqrt{2\pi /P_y}$ and $k$ is the permeability constant.

The radiation term is given as

$$\frac{\partial q_r}{\partial y}=-\frac{4\sigma }{3K}\frac{\partial T^4}{\partial y}\cong -\frac{16\sigma T_s^3}{3K}\frac{{\partial }^{2}T}{\partial y^2}\left(\text{using Rosseland}’\mathrm{s}\text{approximation}\right)$$E36

The appropriate boundary conditions are given as

$$u=0,v=v_w=\frac{\mathit{dh}}{\mathit{dt}},T=T_H\mathit{at}y=h\left(t\right)=H\sqrt{1-\mathit{\alpha t}},$$E37

$$\frac{\partial u}{\partial y}=0,\frac{\partial T}{\partial y}=0,v=0,\mathit{at}y=0,$$E38

On introducing the following dimensionless and similarity variables

$$\begin{array}{l}\begin{array}{l}u=\frac{\alpha H}{2\sqrt{1-\alpha t}}{f}^{\prime }\left(\eta \text{,}t\right),v=-\frac{\alpha H}{2\sqrt{1-\alpha t}}f\left(\eta \text{,}t\right),\eta =\frac{y}{H\sqrt{1-\alpha t}},\hfill \end{array}\hfill \\ \begin{array}{l}\theta =\frac{T-T_0}{T_H-T_0},Ec=\frac{1}{C_p}{\left(\frac{\alpha H}{2\left(1-\alpha t\right)}\right)}^2\hfill \end{array}\hfill \\ \begin{array}{l}\mathrm{Re}=-SA{\left(1-\varphi \right)}^{2.5}=\frac{{\rho }_{nf}H{V}_w}{{\mu }_{nf}},S=\frac{\alpha H^2}{2v_f},Da=\frac{K_p}{H^2},\hfill \end{array}\hfill \\ A_1=\left(1-\varphi \right)+\varphi \frac{{\rho }_s}{{\rho }_f},\mathrm{Pr}=\frac{\mu C_p}{k},\delta =\frac{H}{x},\hfill \\ B_1=\left[\frac{\left({\sigma }_s+\left(m-1\right){\sigma }_f\right)+\left(m-1\right)\left({\sigma }_s-{\sigma }_f\right)\varphi }{\left({\sigma }_s+\left(m-1\right){\sigma }_f\right)-\left(m-1\right)\left({\sigma }_s-{\sigma }_f\right)\varphi }\right],\hfill \\ {\text{A}}_2=\left(1-\varphi \right)+\varphi \frac{{\left(\rho C_p\right)}_s}{{\left(\rho C_p\right)}_f},A_3=\frac{k_{nf}}{k_f},R=\frac{4\sigma T_{\infty }^3}{3kK}\hfill \end{array}$$E39

One arrives at the dimensionless equations

$$\left(1+\frac{1}{\beta }\right)f^{\mathit{iv}}-{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}\left(\eta f^{‴}+3f^{″}+f{f}^{‴}-f^{\prime }{f}^{″}\right)-{\mathit{Ha}}^2{f}^{″}-\frac{1}{\mathit{Da}}{f}^{″}=0$$E40

$$\left(1+\frac{4}{3}R\right){\theta }^{″}+\mathit{PrS}\left(\frac{A_2}{A_3}\right)\left({\theta }^{\prime }f-\eta {\theta }^{\prime }\right)+\frac{Pr\mathit{Ec}}{A_3{\left(1-\varphi \right)}^{2.5}}\left({\left(f^{″}\right)}^2+4{\delta }^2{\left(f^{\prime }\right)}^2\right)=0$$E41

with the boundary conditions as follows

$$f=0,f^{″}=0,{\theta }^{\prime }=0,\mathit{\text{when}}\eta =0,$$E42

$$f=1,f^{\prime }=0,\theta =1,\mathit{\text{when}}\eta =1,$$E43

where m in the above Hamilton Crosser’s model in Eq. (35).

5. Method of solution by homotopy perturbation method

The comparative advantages and the provision of acceptable analytical results with convenient convergence and stability coupled with total analytic procedures of homotopy perturbation method compel us to consider the method for solving the system of nonlinear differential equations in Eqs. (40) and (41) with the boundary conditions in Eq. (42).

5.1 The basic idea of homotopy perturbation method

In order to establish the basic idea behind homotopy perturbation method, consider a system of nonlinear differential equations given as

$$A\left(U\right)-f\left(r\right)=0,r\in \mathrm{\Omega },$$E44

with the boundary conditions

$$B\left(u,\frac{\partial u}{\partial \eta }\right)=0,r\in \mathrm{\Gamma },$$E45

where A is a general differential operator, B is a boundary operator, $f\left(r\right)$ a known analytical function and $\mathrm{\Gamma }$ is the boundary of the domain $\mathrm{\Omega }$.

The operator A can be divided into two parts, which are L and N, where L is a linear operator, N is a non-linear operator. Eq. (44) can be therefore rewritten as follows

$$L\left(u\right)+N\left(u\right)-f\left(r\right)=0.$$E46

By the homotopy technique, a homotopy $U\left(r,p\right):\mathrm{\Omega }\times \left[0,1\right]\to R$ can be constructed, which satisfies

$$H\left(U,p\right)=\left(1-p\right)\left[L\left(U\right)-L\left(U_o\right)\right]+p\left[A\left(U\right)-f\left(r\right)\right]=0,p\in \left[0,1\right],$$E47

or

$$H\left(U,p\right)=L\left(U\right)-L\left(U_o\right)+\mathit{pL}\left(U_o\right)+p\left[N\left(U\right)-f\left(r\right)\right]=0.$$E48

In the above Eqs. (47) and (48), $p\in \left[0,1\right]$ is an embedding parameter, $u_o$ is an initial approximation of equation of Eq. (44), which satisfies the boundary conditions.

Also, from Eq. (47) and Eq. (48), one has

$$H\left(U,0\right)=L\left(U\right)-L\left(U_o\right)=0,$$E49

or

$$H\left(U,0\right)=A\left(U\right)-f\left(r\right)=0.$$E50

The changing process of p from zero to unity is just that of $U\left(r,p\right)$ from $u_o\left(r\right)$ to $u\left(r\right)$. This is referred to homotopy in topology. Using the embedding parameter p as a small parameter, the solution of Eqs. (47) and Eq. (48) can be assumed to be written as a power series in p as given in Eq. (51)

$$U=U_o+{\mathit{pU}}_1+p^2{U}_2+\dots$$E51

It should be pointed out that of all the values of p between 0 and 1, p = 1 produces the best result. Therefore, setting $p=1$, results in the approximation solution of Eq. (42)

$$u=\underset{p\to 1}{lim}U=U_o+U_1+U_2+\dots$$E52

The basic idea expressed above is a combination of homotopy and perturbation method. Hence, the method is called homotopy perturbation method (HPM), which has eliminated the limitations of the traditional perturbation methods. On the other hand, this technique can have full advantages of the traditional perturbation techniques. The series Eq. (29) is convergent for most cases.

5.2 Application of the homotopy perturbation method to the fluid flow problem

According to homotopy perturbation method (HPM), one can construct an homotopy for Eq. (36)–(39) as

$$H_1\left(p,\eta \right)=\left(1-p\right)\left[\left(1+\frac{1}{\beta }\right)f^{\mathit{iv}}\right]+p\left[\begin{array}{l}\left(1+\frac{1}{\beta }\right)f^{\mathit{iv}}-{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}\left(\begin{array}{l}\eta f^{‴}+3f^{″}\\ +f{f}^{‴}-f^{\prime }{f}^{″}\end{array}\right)\\ -{\mathit{Ha}}^2{f}^{″}-\frac{1}{\mathit{Da}}{f}^{″}\end{array}\right]=0,$$E53

$$H_2\left(p,\eta \right)=\left(1-p\right)\left[\left(1+\frac{4}{3}R\right){\theta }^{″}\right]+p\left[\begin{array}{l}\left(1+\frac{4}{3}R\right){\theta }^{″}+\mathit{PrS}\left(\frac{A_2}{A_3}\right)\left({\theta }^{\prime }f-\eta {\theta }^{\prime }\right)\\ +\frac{Pr\mathit{Ec}}{A_3{\left(1-\varphi \right)}^{2.5}}\left({\left(f^{″}\right)}^2+4{\delta }^2{\left(f^{\prime }\right)}^2\right)\end{array}\right]=0,$$E54

Taking power series of velocity, temperature and concentration fields, gives

$$f=f_0+{\mathit{pf}}_1+p^2{f}_2+p^3{f}_3+\dots$$E55

and

$$\theta ={\theta }_0+p{\theta }_1+p^2{\theta }_2+p^3{\theta }_3+\dots$$E56

Substituting Eqs. (55) and (56) into Eq. (53) and (54) as well as the boundary conditions in Eq. (42), and grouping like terms based on the power of p, the fluid flow velocity equation is given as:

Zeroth-order equations

$$p^0:f_0^{\mathit{iv}}\left(\eta \right)+\frac{1}{\beta }{f}_0^{\mathit{iv}}\left(\eta \right)=0,$$E57

$$p^0:\left(1+\frac{4}{3}R\right){\theta }_0^{″}=0,$$E58

First-order equations

$$\begin{array}{l}{p}^1:\frac{1}{\beta }{f}_1^{\mathit{iv}}\left(\eta \right)+f_1^{\mathit{iv}}\left(\eta \right)-{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}\eta f_0\left(\eta \right)-\frac{1}{\mathit{Da}}{f}_0^{″}\left(\eta \right)-{\mathit{Ha}}^2{f}_0^{″}\left(\eta \right)\\ -3{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{f}_0^{″}\left(\eta \right)-{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{f}_0^{\prime }\left(\eta \right)f_0^{″}\left(\eta \right)+{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{f}_0\left(\eta \right)f_0^{‴}\left(\eta \right)=0,\end{array}$$E59

$$p^1:\left(1+\frac{4}{3}R\right){\theta }_1^{″}+\mathit{PrS}\left(\frac{A_2}{A_3}\right)\left({\theta }_0^{\prime }{f}_0-\eta {\theta }_0^{\prime }\right)+\frac{Pr\mathit{Ec}}{A_3{\left(1-\varphi \right)}^{2.5}}\left({\left(f_0^{″}\right)}^2+4{\delta }^2{\left(f_0^{\prime }\right)}^2\right)=0$$E60

Second-order equations

$$\begin{array}{l}{p}^2:\frac{1}{\beta }{f}_2^{\mathit{iv}}\left(\eta \right)+f_2^{\mathit{iv}}\left(\eta \right)-{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}\eta f_1\left(\eta \right)-\frac{1}{\mathit{Da}}{f}_1^{″}\left(\eta \right)-{\mathit{Ha}}^2{f}_1^{″}\left(\eta \right)-3{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{f_2}^{\prime \prime }\left(\eta \right)\\ -{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{f}_1^{\prime }\left(\eta \right)f_0^{″}\left(\eta \right)-{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{f}_0^{\prime }\left(\eta \right)f_1^{″}\left(\eta \right)+{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{f}_1\left(\eta \right)f_0^{‴}\left(\eta \right)\\ +{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{f}_0\left(\eta \right)f_1^{‴}\left(\eta \right)=0,\end{array}$$E61

$$p^2:\left(1+\frac{4}{3}R\right){\theta }_2^{″}+\mathit{PrS}\left(\frac{A_2}{A_3}\right)\left({\theta }_1^{\prime }{f}_0+{\theta }_0^{\prime }{f}_1-\eta {\theta }_1^{\prime }\right)+\frac{2Pr\mathit{Ec}}{A_3{\left(1-\varphi \right)}^{2.5}}\left(f_0^{″}{f}_1^{″}+4{\delta }^2{f_0}^{\prime }{f}_1^{\prime }\right)=0$$E62

the boundary conditions are

$$\begin{array}{l}{f}_0=f_1=f_2=0,{f_0}^{\prime \prime }={f_1}^{\prime \prime }={f_2}^{\prime \prime }=0,{{\theta }_0}^{\prime }={{\theta }_1}^{\prime }={{\theta }_2}^{\prime }=0,\mathit{\text{when}}\eta =0,\\ f_0=1,f_1=f_2=0,{f_0}^{\prime }={f_1}^{\prime }={f_2}^{\prime }=0,{\theta }_0=1,{\theta }_1={\theta }_2=0,\mathit{\text{when}}\eta =1,\end{array}$$E63

In a similar way, the higher orders problems are obtained.

On solving Eqs. (57), (61) and (64) with their corresponding boundary conditions, we arrived at

$$f_0\left(\eta \right)=\frac{1}{2}\left(3\eta -{\eta }^3\right)$$E64

$$f_1\left(\eta \right)=-\frac{1}{6720\left(1+\beta \right)}\left(\begin{array}{l}\left(168\left(\frac{1}{\mathit{Da}}\right)\beta +168{\mathit{Ha}}^2\beta +419{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}\beta \right)\eta \\ -\left(336\left(\frac{1}{\mathit{Da}}\right)\beta +336{\mathit{Ha}}^2\beta +873{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}\beta \right){\eta }^3\\ +\left(168\left(\frac{1}{\mathit{Da}}\right)\beta +168{\mathit{Ha}}^2\beta +504{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}\beta \right){\eta }^5\\ -28{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{\mathit{\beta \eta }}^6-24{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{\mathit{\beta \eta }}^7\\ +2{\mathit{SA}}_1{\left(1-\varphi \right)}^{2.5}{\mathit{\beta \eta }}^8\end{array}\right)$$E65

E66

In the same manner, the energy equations are solved. Following the definition of the homotopy perturbation method as presented in Eq. (52), one could write the solution of the fluid flow equation as

E67

6. Example 3: homotopy perturbation method to dynamic behavior of piezoelectric nanobeam embedded in linear and nonlinear elastic Foundation in a thermal-magnetic environment

Consider a nanobeam embedded in linear and nonlinear elastic media as shown in Figure 3. The nanobeam is subjected to stretching effects and resting on Winkler, Pasternak and nonlinear elastic media in a thermo-magnetic environment as depicted in the figure.

Figure 3.

A nanobeam embedded in linear and nonlinear elastic media (note: Only the bottom side of the elastic media is shown).

Following the nonlocal theory and Euler-Bernoulli theorem, the governing equation of the structure is developed as

$$\begin{array}{l}\mathit{EI}\left(\frac{{\partial }^4\overline{w}}{\partial {\overline{x}}^4}\right)+\rho A_c\frac{{\partial }^2}{\partial {\overline{t}}^2}\left[\overline{w}-{\left(e_{0}a\right)}^2\frac{{\partial }^2\overline{w}}{\partial {\overline{x}}^2}\right]+k_w\left[\overline{w}-{\left(e_{0}a\right)}^2\frac{{\partial }^2\overline{w}}{\partial {\overline{x}}^2}\right]-k_p\frac{{\partial }^2}{\partial {\overline{x}}^2}\left[\overline{w}-{\left(e_{0}a\right)}^2\frac{{\partial }^2\overline{w}}{\partial {\overline{x}}^2}\right]\\ +k_2\left[{\overline{w}}^2-{\left(e_{0}a\right)}^2\frac{{\partial }^2\left({\overline{w}}^2\right)}{\partial {\overline{x}}^2}\right]+k_3\left[{\overline{w}}^3-{\left(e_{0}a\right)}^2\frac{{\partial }^2\left({\overline{w}}^3\right)}{\partial {\overline{x}}^2}\right]-\eta A_c{H}_{\overline{x}}^2\frac{{\partial }^2}{\partial {\overline{x}}^2}\left[\overline{w}-{\left(e_{0}a\right)}^2\frac{{\partial }^2\overline{w}}{\partial {\overline{x}}^2}\right]\\ +\left({\mathit{EA}}_c\frac{{\alpha }_{\overline{x}}\mathrm{\Delta }T}{1-2\nu }\right)\frac{{\partial }^2}{\partial {\overline{x}}^2}\left[\overline{w}-{\left(e_{0}a\right)}^2\frac{{\partial }^2\overline{w}}{\partial {\overline{x}}^2}\right]-\left[\left(\frac{{\mathit{EA}}_c}{2L}{\int }_0^L{\left(\frac{\partial \overline{w}}{\partial \overline{x}}\right)}^{2}d\overline{x}\right)\left(\frac{{\partial }^2\overline{w}}{\partial {\overline{x}}^2}-{\left(e_{0}a\right)}^2\frac{{\partial }^4\overline{w}}{\partial {\overline{x}}^4}\right)\right]=0\end{array}$$E68

It is assumed that the midpoint of the nanobeam is subjected to the following initial conditions

$$\overline{w}\left(\overline{x},0\right)={\overline{w}}_o,\frac{\partial \overline{w}\left(\overline{x},0\right)}{\partial \overline{t}}=0$$E69

The following boundary conditions for the multi-walled nanotubes for simply supported nanotube is given,

$$\overline{w}\left(0,\overline{t}\right)=0,\frac{{\partial }^2\overline{w}\left(0,\overline{t}\right)}{{\partial }^2\overline{x}}=0,\overline{w}\left(L,\overline{t}\right)=0,\frac{{\partial }^2\overline{w}\left(L,\overline{t}\right)}{{\partial }^2\overline{x}}=0.$$E70

Using the following adimensional constants and variables

$$\begin{array}{l}x=\frac{\overline{x}}{L};w=\frac{\overline{w}}{r};t=\sqrt{\frac{\mathit{EI}}{\rho A_c{L}^4}};r=\sqrt{\frac{I}{A_c}};h=\frac{e_{0}a}{L};{\alpha }_t^d=\frac{N_{\mathit{\text{thermal}}}{L}^2}{\mathit{EI}};A=\frac{{\overline{w}}_o}{r}\\ K_w=\frac{k_w{L}^4}{\mathit{EI}};K_p=\frac{k_{p}L2}{\mathit{EI}};{\mathit{Ha}}_m=\frac{\eta A_c{H}_{\overline{x}}^2{L}^2}{\mathit{EI}};K_2^d=\frac{k_2{\mathit{rL}}^4}{\mathit{EI}};K_3^d=\frac{k_3{r}^2{L}^4}{\mathit{EI}}.\end{array}$$E71

The adimensional form of the governing equation of motion for the nanobeam is given as

$$\begin{array}{l}\left[1+K_p{h}^2+H{a}_m{h}^2-{\alpha }_t^d{h}^2+\frac{h^2}{2}\underset{0}{\overset{1}{{\displaystyle \int }}}{\left(\frac{\partial w}{\partial x}\right)}^{2}dx\right]\hfill \\ \frac{{\partial }^{4}w}{\partial x^4}+\left[{\alpha }_t^d-K_w{h}^2-K_p-H{a}_m-\frac{1}{2}\underset{0}{\overset{1}{{\displaystyle \int }}}{\left(\frac{\partial w}{\partial x}\right)}^{2}dx\right]\frac{{\partial }^{2}w}{\partial x^2}\hfill \\ +K_{w}w+\frac{{\partial }^{2}w}{\partial t^2}-h^2\frac{{\partial }^{4}w}{\partial x^2\partial t^2}+K_2^d\left[w^2-h^2\frac{{\partial }^2\left(w^2\right)}{\partial x^2}\right]+K_3^d\left[w^3-h^2\frac{{\partial }^2\left(w^3\right)}{\partial x^2}\right]=0\hfill \end{array}$$E72

And the boundary conditions become

$$w\left(0,t\right)=0,\frac{{\partial }^{2}w\left(0,t\right)}{{\partial }^{2}x}=0,w\left(1,t\right)=0,\frac{{\partial }^{2}w\left(1,t\right)}{{\partial }^{2}x}=0.$$E73

6.1 Solution methodology: Galerkin decomposition and homotopy perturbation methods

The method of solution for the governing equation includes Galerkin decomposition and homotopy perturbation methods. As the name implies the Galerkin decomposition method is used to decompose the governing partial differential equation of motion can be separated into spatial and temporal parts. The resulting temporal equations are solved using homotopy perturbation method.

The procedures for the analysis of the equations are given in the proceeding sections as follows:

6.1.1 Galerkin decomposition method

With the application of Galerkin decomposition procedure, the governing partial differential equations of motion can be separated into spatial and temporal parts of the lateral displacement function as

$$w\left(x,t\right)=\varphi \left(x\right)q\left(t\right)$$E74

Using one-parameter Galerkin decomposition procedure, one arrives at

$$\underset{0}{\overset{1}{\int }}R\left(x,t\right)\varphi \left(x\right)\mathit{dx}=0$$E75

where $R\left(x,t\right)$ is the governing equation of motion for nanobeam i.e.

$$\begin{array}{l}R\left(x\text{,}t\right)=\left[1+K_p{h}^2+H{a}_m{h}^2-{\alpha }_t^d{h}^2+\frac{h^2}{2}\underset{0}{\overset{1}{{\displaystyle \int }}}{\left(\frac{\partial w}{\partial x}\right)}^{2}dx\right]\frac{{\partial }^{4}w}{\partial x^4}\hfill \\ +\left[{\alpha }_t^d-K_w{h}^2-K_p-H{a}_m-\frac{1}{2}\underset{0}{\overset{1}{{\displaystyle \int }}}{\left(\frac{\partial w}{\partial x}\right)}^{2}dx\right]\frac{{\partial }^{2}w}{\partial x^2}\hfill \\ +K_{w}w+\frac{{\partial }^{2}w}{\partial t^2}-h^2\frac{{\partial }^{4}w}{\partial x^2\partial t^2}+K_2^d\left[w^2-h^2\frac{{\partial }^2\left(w^2\right)}{\partial x^2}\right]+K_3^d\left[w^3-h^2\frac{{\partial }^2\left(w^3\right)}{\partial x^2}\right]=0\hfill \end{array}$$E76

where $\varphi \left(x\right)$ is the basis or trial or comparison function or normal function, which must satisfy the boundary conditions in Eq. (73), and $q\left(t\right)$ is the temporal part (time-dependent function).

Substituting Eqs. (75) into (74), then multiplying both sides of the resulting equation by $\varphi \left(x\right)$ and integrating it for the domain of (0,1), we have

$$\frac{d^{2}q\left(t\right)}{{\mathit{dt}}^2}+{\lambda }_{1}q\left(t\right)+{\lambda }_2{q}^2\left(t\right)+{\lambda }_3{q}^3\left(t\right)=0$$E77

where

$${\lambda }_1=\frac{{\overline{\lambda }}_1}{{\overline{\lambda }}_0};{\lambda }_2=\frac{{\overline{\lambda }}_2}{{\overline{\lambda }}_0};{\lambda }_3=\frac{{\overline{\lambda }}_3}{{\overline{\lambda }}_0};$$E78

$${\overline{\lambda }}_0=\underset{0}{\overset{1}{\int }}\left({\varphi }^2-h^2\varphi \frac{{\partial }^2\varphi }{\partial x^2}\right)\mathit{dx}$$E79

$${\overline{\lambda }}_1=\underset{0}{\overset{1}{\int }}\left(K_w{\varphi }^2+\left(1+K_p{h}^2+{\mathit{Ha}}_m{h}^2-{\alpha }_t^d{h}^2\right)\varphi \frac{{\partial }^4\varphi }{\partial x^4}+\left({\alpha }_t^d-K_w{h}^2-K_p-{\mathit{Ha}}_m\right)\varphi \frac{{\partial }^2\varphi }{\partial x^2}\right)\mathit{dx}$$E80

$${\overline{\lambda }}_2=\underset{0}{\overset{1}{\int }}{K}_2^d\left({\varphi }^3-h^2\varphi \frac{{\partial }^2\left({\varphi }^2\right)}{\partial x^2}\right)\mathit{dx}$$E81

$${\overline{\lambda }}_3=\underset{0}{\overset{1}{\int }}{K}_3^d\left({\varphi }^4-h^2\varphi \frac{{\partial }^2\left({\varphi }^4\right)}{\partial x^2}\right)\mathit{dx}+\frac{h^2}{2}\underset{0}{\overset{1}{\int }}{\left(\frac{\partial \varphi }{\partial x}\right)}^2\mathit{dx}\underset{0}{\overset{1}{\int }}\varphi \frac{{\partial }^2\varphi }{\partial x^2}\mathit{dx}-\frac{1}{2}\underset{0}{\overset{1}{\int }}{\left(\frac{\partial \varphi }{\partial x}\right)}^2\mathit{dx}\underset{0}{\overset{1}{\int }}\varphi \frac{{\partial }^4\varphi }{\partial x^4}\mathit{dx}$$E82

The initial conditions are given as

$$q\left(0\right)=A,\frac{\mathit{dq}\left(0\right)}{\mathit{dt}}=0$$E83

A is the maximum vibration amplitude of the structure.

From the initial conditions in Eq. (83), one can write the initial approximation, $u_o$ as

$$u_o=\mathit{\text{Acos}}\left(\mathit{\omega t}\right)$$E84

Eq. (22) satisfies the initial conditions in Eq. (83).

The homotopy perturbation representation of Eq. (77) is

$$H\left(q,p\right)=\left[\frac{d^{2}q}{{\mathit{dt}}^2}+{\lambda }_{1}q\right]-\left[\frac{d^2{u}_o}{{\mathit{dt}}^2}+{\lambda }_1{u}_o\right]+p\left[\frac{d^2{u}_o}{{\mathit{dt}}^2}+{\lambda }_1{u}_o\right]+p\left({\lambda }_2{q}^2+{\lambda }_3{q}^3\right)=0$$E85

From the procedure of homotopy perturbation method, assuming that the solution of Eq. (77) takes the form of:

$$q=q_0+{\mathit{pq}}_1+p^2{q}_2+p^3{q}_3+\dots ,$$E86

On substituting Eqs. (86) into the homotopy Eq. (85)

$$\begin{array}{l}H\left(q,p\right)=\left[\frac{d^2\left(q_0+{\mathit{pq}}_1+p^2{q}_2+p^3{q}_3+\dots \right)}{{\mathit{dt}}^2}+{\lambda }_1\left(q_0+{\mathit{pq}}_1+p^2{q}_2+p^3{q}_3+\dots \right)\right]\\ -\left[\frac{d^2{u}_o}{{\mathit{dt}}^2}+{\lambda }_1{u}_0\right]+p\left[\frac{d^2{u}_o}{{\mathit{dt}}^2}+{\lambda }_1{u}_0\right]+p\left(\begin{array}{l}{\lambda }_2{\left(q_0+{\mathit{pq}}_1+p^2{q}_2+p^3{q}_3+\dots \right)}^2\\ +{\lambda }_3{\left(q_0+{\mathit{pq}}_1+p^2{q}_2+p^3{q}_3+\dots \right)}^3\end{array}\right)=0\end{array}$$E87

rearranging the coefficients of the terms with identical powers of p, one obtains series of linear differential equations as.

Zero-order equation

$$p^0:\left[\frac{d^2{q}_0}{{\mathit{dt}}^2}+{\lambda }_1{q}_0\right]-\left[\frac{d^2{u}_o}{{\mathit{dt}}^2}+{\lambda }_1{u}_o\right]=0$$E88

with the conditions

$$q_0\left(0\right)=A\text{and}\frac{{\mathit{dq}}_0\left(0\right)}{\mathit{dt}}=0$$E89

First-order equation

$$p^1:\frac{d^2{q}_1}{{\mathit{dt}}^2}+{\lambda }_1{q}_0+\frac{d^2{u}_o}{{\mathit{dt}}^2}+{\lambda }_1{u}_o+{\lambda }_2{q}_0^2+{\lambda }_3{q}_0^3=0$$E90

with corresponding initial conditions

$$q_1\left(0\right)=0\text{and}\frac{{\mathit{dq}}_1\left(0\right)}{\mathit{dt}}=0$$E91

Second-order equation

$$p^2:\frac{d^2{q}_2}{{\mathit{dt}}^2}+{\lambda }_1{q}_2+2{\lambda }_2{q}_0{q}_1+3{\lambda }_3{q}_0^2{q}_1=0$$E92

with corresponding initial conditions

$$q_2\left(0\right)=0\text{and}\frac{{\mathit{dq}}_2\left(0\right)}{\mathit{dt}}=0$$E93

The solution of the zero-order is given by.

From Eq. (27), we have

$$q_0=\mathit{\text{Acos}}\left(\mathit{\omega t}\right)$$E94

On substituting Eq. (94) into Eq. (90) and using trigonometric identities, after the colllection of like terms, one arrives at

$$\frac{d^2{q}_1}{{\mathit{dt}}^2}+{\lambda }_1{q}_1+A\left({\lambda }_1-{\omega }^2+\frac{3}{4}{A}^2\lambda \right)\mathit{cos}\left(\mathit{\omega t}\right)+\frac{A^2{\lambda }_2}{2}\mathit{cos}\left(2\mathit{\omega t}\right)+\frac{A^3{\lambda }_3}{4}\mathit{cos}\left(3\mathit{\omega t}\right)+\frac{A^2{\lambda }_2}{2}=0$$E95

The solution of the above Eq. (95) provides

$$\begin{array}{l}{q}_1\left(t\right)=\left[\begin{array}{l}A\left({\lambda }_1-{\omega }^2+\frac{3}{4}{A}^2\lambda \right)\left(\frac{{\lambda }_1}{{\omega }^2-{\lambda }_1^2}\right)\cos \left(\omega t\right)+\frac{A^2{\lambda }_2}{2}\left(\frac{{\lambda }_1}{4{\omega }^2-{\lambda }_1^2}\right)\cos \left(2\omega t\right)\hfill \\ +\frac{A^3{\lambda }_3}{4}\left(\frac{{\lambda }_1}{9{\omega }^2-{\lambda }_1^2}\right)\cos \left(3\omega t\right)+\frac{A^2{\lambda }_2}{2}\hfill \end{array}\right]\hfill \\ +[A\left({\lambda }_1-{\omega }^2+\frac{3}{4}{A}^2\lambda \right)\left(\frac{{\lambda }_1}{{\lambda }_1^2-{\omega }^2}\right)+\frac{A^2{\lambda }_2}{2}\left(\frac{{\lambda }_1}{{\lambda }_1^2-4{\omega }^2}\right)\hfill \\ +\frac{A^3{\lambda }_3}{4}\left(\frac{{\lambda }_1}{{\lambda }_1^2-9{\omega }^2}\right)+\frac{A^2{\lambda }_2}{2{\lambda }_1}]\cos \left(\alpha t\right)\hfill \end{array}$$E96

Based on the procedure of HPM, setting $p=1$,

$$q\left(t\right)=\underset{p\to 1}{lim}q\left(t\right)=\underset{p\to 1}{lim}\left[q_0+{\mathit{pq}}_1+p^2{q}_2+p^3{q}_3+\dots \right]=q_0+q_1+q_2+q_3+\dots$$E97

On substituting Eqs. (94) and (96) into Eq. (97), the result is

$$\begin{array}{l}q\left(t\right)=Acos\left(\omega t\right)+\left[\begin{array}{l}A\left({\lambda }_1-{\omega }^2+\frac{3}{4}{A}^2\lambda \right)\left(\frac{{\lambda }_1}{{\omega }^2-{\lambda }_1^2}\right)\cos \left(\omega t\right)+\frac{A^2{\lambda }_2}{2}\left(\frac{{\lambda }_1}{4{\omega }^2-{\lambda }_1^2}\right)\cos \left(2\omega t\right)\hfill \\ +\frac{A^3{\lambda }_3}{4}\left(\frac{{\lambda }_1}{9{\omega }^2-{\lambda }_1^2}\right)\cos \left(3\omega t\right)+\frac{A^2{\lambda }_2}{2}\hfill \end{array}\right]\hfill \\ +[A\left({\lambda }_1-{\omega }^2+\frac{3}{4}{A}^2\lambda \right)\left(\frac{{\lambda }_1}{{\lambda }_1^2-{\omega }^2}\right)+\frac{A^2{\lambda }_2}{2}\left(\frac{{\lambda }_1}{{\lambda }_1^2-4{\omega }^2}\right)\hfill \\ +\frac{A^3{\lambda }_3}{4}\left(\frac{{\lambda }_1}{{\lambda }_1^2-9{\omega }^2}\right)+\frac{A^2{\lambda }_2}{2{\lambda }_1}]\cos \left({\lambda }_{1}t\right)+\dots \hfill \end{array}$$E98

In order to find the natural frequency, $\omega$, the secular term must be eliminated. In order to do this, set the coefficient of $\mathit{cos}\left({\lambda }_{1}t\right)$ to zero.

$$A\left({\lambda }_1-{\omega }^2+\frac{3}{4}{A}^2\lambda \right)\left(\frac{{\lambda }_1}{{\lambda }_1^2-{\omega }^2}\right)+\frac{A^2{\lambda }_2}{2}\left(\frac{{\lambda }_1}{{\lambda }_1^2-4{\omega }^2}\right)+\frac{A^3{\lambda }_3}{4}\left(\frac{{\lambda }_1}{{\lambda }_1^2-9{\omega }^2}\right)+\frac{A^2{\lambda }_2}{2{\lambda }_1}=0$$E99

After simplification of Eq. (99), we have

$$\begin{array}{l}\left(\frac{A{\lambda }_2}{2{\lambda }_1^2}-1\right){\omega }^6+A\left[{\lambda }_1^2\left(13-\frac{49A{\lambda }_2}{2}\right)-36{\lambda }_1+\frac{9A{\lambda }_2}{2}-26{\lambda }_{3}A\right]{\omega }^4\\ A\left[{\lambda }_1^4+13{\lambda }_1^3-\left(2A{\lambda }_2-11{\lambda }_3{A}^2\right){\lambda }_1^2\right]{\omega }^2+{\lambda }_1^{4}A\left({\lambda }_1+{\lambda }_3{A}^2\right)=0\end{array}$$E100

The sextic equation can be written as

$$\begin{array}{l}\left(\frac{A{\lambda }_2}{2{\lambda }_1^2}-1\right){\omega }^6+A\left[{\lambda }_1^2\left(13-\frac{49A{\lambda }_2}{2}\right)-36{\lambda }_1+\frac{9A{\lambda }_2}{2}-26{\lambda }_{3}A\right]{\omega }^4\\ A\left[{\lambda }_1^4+13{\lambda }_1^3-\left(2A{\lambda }_2-11{\lambda }_3{A}^2\right){\lambda }_1^2\right]{\omega }^2+{\lambda }_1^{4}A\left({\lambda }_1+{\lambda }_3{A}^2\right)=0\end{array}$$E101

Eq. (101) can be written as

$${\chi }_1{\omega }^6+{\chi }_2{\omega }^4+{\chi }_3{\omega }^2+{\chi }_4=0$$E102

where

$${\chi }_1=\left(\frac{A{\lambda }_2}{2{\lambda }_1^2}-1\right),{\chi }_2=A\left[{\lambda }_1^2\left(13-\frac{49A{\lambda }_2}{2}\right)-36{\lambda }_1+\frac{9A{\lambda }_2}{2}-26{\lambda }_{3}A\right]$$

$${\chi }_3=A\left[{\lambda }_1^4+13{\lambda }_1^3-\left(2A{\lambda }_2-11{\lambda }_3{A}^2\right){\lambda }_1^2\right],{\chi }_4={\lambda }_1^{4}A\left({\lambda }_1+{\lambda }_3{A}^2\right)=0$$

The roots of the sextic equation are

$$\begin{array}{l}{\omega }_1=\sqrt{\begin{array}{l}\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)+\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\\ +\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)-\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}-\frac{{\chi }_2}{3{\chi }_1}\end{array}}\\ \end{array}$$E103

$$\begin{array}{l}{\omega }_2=-\sqrt{\begin{array}{l}\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)+\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\\ +\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)-\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}-\frac{{\chi }_2}{3{\chi }_1}\end{array}}\\ \end{array}$$E104

$$\begin{array}{l}{\omega }_3=\sqrt{\begin{array}{l}\frac{-1}{2{\chi }_1}\left[\begin{array}{l}\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)+\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\\ +\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)-\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\end{array}\right]\\ +\frac{\sqrt{-3}}{2{\chi }_1}\left[\begin{array}{l}\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)+\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\\ -\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)-\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\end{array}\right]-\frac{{\chi }_2}{3{\chi }_1}\end{array}}\\ \end{array}$$E105

$$\begin{array}{l}{\omega }_4=-\sqrt{\begin{array}{l}\frac{-1}{2{\chi }_1}\left[\begin{array}{l}\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)+\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\\ +\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)-\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\end{array}\right]\\ +\frac{\sqrt{-3}}{2{\chi }_1}\left[\begin{array}{l}\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)+\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\\ -\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)-\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\end{array}\right]-\frac{{\chi }_2}{3{\chi }_1}\end{array}}\\ \end{array}$$E106

$$\begin{array}{l}{\omega }_5=\sqrt{\begin{array}{l}\frac{-1}{2{\chi }_1}\left[\begin{array}{l}\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)+\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\\ +\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)-\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\end{array}\right]\\ \frac{-\sqrt{-3}}{2{\chi }_1}\left[\begin{array}{l}\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)+\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\\ -\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)-\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\end{array}\right]-\frac{{\chi }_2}{3{\chi }_1}\end{array}}\\ \end{array}$$E107

$${\omega }_6=-\sqrt{\begin{array}{l}\frac{-1}{2{\chi }_1}\left[\begin{array}{l}\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)+\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\\ +\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)-\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\end{array}\right]\\ \frac{-\sqrt{-3}}{2{\chi }_1}\left[\begin{array}{l}\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)+\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\\ -\sqrt[3]{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)-\sqrt{{\left(\frac{{\chi }_3}{3{\chi }_1}-\frac{{\lambda }_2^2}{9{\chi }_1^2}\right)}^3+{\left(\frac{-{\chi }_2^3}{27{\chi }_1^3}+\frac{{\chi }_2{\chi }_3}{6{\chi }_1^2}-\frac{{\chi }_4}{2{\chi }_1}\right)}^2}}\end{array}\right]-\frac{{\chi }_2}{3{\chi }_1}\end{array}}$$E108

7. Conclusion

In this chapter, the applications of regular and homotopy perturbation methods to thermal, fluid flow and dynamic behaviors of engineering systems have been presented. Regular perturbation was used in the first example to developed approximate analytical solutions for thermal behavior of convective-radiative fin with end cooling and thermal contact resistance. In the second example, homotopy perturbation method utilized to study squeezing flow and heat transfer of Casson nanofluid between two parallel plates embedded in a porous medium under the influences of slip, Lorentz force, viscous dissipation and thermal radiation. The same method was used in the third example to analyze the dynamic behavior of piezoelectric nanobeam embedded in linear and nonlinear elastic foundations operating in a thermal-magnetic environment. It is hoped that the vivid presentation and applications of these perturbation methods in this chapter will advance better understanding of methods especially for real world applications.