# Measuring the Pollutants in a System of Three Interconnecting Lakes by the Semianalytical Method

Authors: Indranil Ghosh; M. S. H. Chowdhury; Suazlan Mt Aznam; M. M. Rashid
Journal: Journal of Applied Mathematics
 (2021)

---

## Abstract

Pollution has become an intense danger to our environment. The lake pollution model is formulated into the three-dimensional system of differential equations with three instances of input. In the present study, the new iterative method (NIM) was applied to the lake pollution model with three cases called impulse input, step input, and sinusoidal input for a longer time span. The main feature of the NIM is that the procedure is very simple, and it does not need to calculate any special type of polynomial or multipliers such as Adomian polynomials and Lagrange’s multipliers. Comparisons with the Adomian decomposition method (ADM) and the well-known purely numerical fourth-order Runge-Kutta method (RK4) suggest that the NIM is a powerful alternative for differential equations providing more realistic series solutions that converge very rapidly in real physical problems.


---

## Body

## 1. Introduction

Over the years, differential equations are employed by researchers that care about what is going on in the surroundings, for example, one can see the following ref. [1–10]. Therefore, differential equations play a prominent role to solve physical problems using numerical or analytical methods. There has recently been much attention devoted to the search for better and more efficient solution methods for determining a solution, approximate or exact, analytical, or numerical, to physical models. Finding exact or approximate solutions of these nonlinear equations is interesting and important. Several semianalytic methods including the Adomian decomposition method (ADM) [11], homotopy perturbation method (HPM) [12], and variational iteration method (VIM) [13] have been used for solving mathematical models of differential equations.Finding accurate and analytical procedures for solving lake pollution model has been a vital research undertaking. The system of lake pollution models has been solved first time by Biazar et al. [14] using a semianalytic technique called the Adomian decomposition method (ADM). So often, it is very hard to ascertain the Adomian polynomials involved with ADM. To overcome the shortcomings of ADM, many analytical techniques have been successfully implemented by various groups of researchers. In particular, the lake pollution system was solved by Ganji et al. [15] via homotopy-perturbation method (HPM). Recently, Ahmed and Khan [16] have investigated dynamics of system of polluted lakes through the numerical procedure fractional derivative approach. Numerous authors used some of the approximate numerical techniques to obtain approximate analytical solutions such as the Bessel matrix method [17], Haar wavelet collocation method [18], Differential Transform method (DTM) [19], Laplace-Padé Differential transform method (LPDTM) [20], and variational iteration method (VIM) [21].Another prominent semianalytical technique which has been demonstrated to be a lot more straightforward and efficient than the abovementioned methods is called the new iterative method (NIM), first proposed by Varsha Daftardar-Gejji and Jafari [22] . Very recently, NIM was employed for solving Fokker-Planck equations [23], chemical kinetic equations [24], and Klein-Gordon equations [25]. Some other applications of the new iterative method (NIM) on nonlinear problems can be found in [26–29]. In this paper, the NIM motivates us to solve the lake pollution systems.


## 2. Mathematical Modeling of the Lake Pollution System

Pollutants have different kinds of characteristic property, and nonpolar pollutants never mix with lake water, whereas polar one makes a reaction with lake water; as a result, mild to severe contamination takes place [30]. Figure 1 represents the three lakes with interconnecting passage with each other. We contemplate every lake to be a big compartment and interconnecting passage as pipes between the compartments. In the figure, each arrow pointed out the flow direction of the pollutant through the channels and qtrepresents the rate of the pollutant that enters the lake per unit time. The volume of water vi in each lake remains constants, and amount of the pollutant is denoted by zit, constantly distributed in each lake at any time t≥0 where i=1,2,3. Here, the concentration of the contaminant in each lake is represented by cit=zit/vi. We consider the constant flow rate is Fji and the motion of the flow is in the direction from lake i to lake j, and there is no reverse flow. The flux rate is symbolized by Rjit and can be defined as Rjit=Fjicit=Fjizit/vi .Figure 1
Interconnecting channels of the system of three lakes [14].Now, the rate of change ofcontaminant=inputrate in every lake is the output rate in every lake. For the volume of every lake to stay steady, the stream rate into every lake must adjust the progression of the lake. So, we acquire the accompanying conditions between (1), (2), and (3):(1)LakeI.F13=F31+F21.(2)LakeII.F32=F21.(3)LakeIII.F31+F32=F13.The following system of differential equations yielded as the lake pollution model:(4)dz1dt=F13v3z3t+qt−F31v1z1t−F21v1z1t,dz2dt=F21v1z1t−F32v2z2t,dz3dt=F31v1z1t+F32v2z2t−F13v3z3t.With initial conditions,(5)z10=z20=z30=0.The parametric values arev1=2900mi3,v2=850mi3,v3=1180mi3,(6)F21=18mi3/year,F32=18mi3/year,F31=20mi3/year,F13=38mi3/year.


## 3. Analysis of the Method

### 3.1. New Iterative Method (NIM)

In this section, we have discussed the new iterative method (NIM) [31, 32] as follows. Let us consider the general equation:(7)uy=Nuy+fy,where N is a nonlinear operator from a Banach space B⟶B, and f is a known function.Lety=y1,y2,⋯..,yn be a solution ofu for equation (7) having the series form.(8)uy=∑i=0∞uiy.The nonlinear operatorN can be written as(9)N∑i=0∞ui=Nu0+∑i=1∞N∑j=0iuj−N∑j=0i−1uj.Then, series solution becomes(10)∑i=1∞ui=f+Nu0+∑i=1∞N∑j=0iuj−N∑j=0i−1uj.We define the consecutive relation:(11)u0=f,(12)u1=Nu0,(13)u2=Nu0+u1−Nu0,(14)um+1=Nu0+⋯+um−Nu0+⋯+um−1,m=1,2,⋯Therefore,(15)u=f+∑i=1∞ui.


### 3.2. Stability Analysis

LetX,Y be Banach spaces and F:X⟶Y a map. LX,Y denotes the set of all linear maps from X to Y. LX,Y is also a Banach space.Definition 1 (see [33, 34]).
F is said to be Fréchet differentiable at x∈X if there exists a continuous linear map A:X⟶Y such that(16)Fx+h−Fx=Ah+wx,h,where limh⟶0wx,h/h=0.A is called the Fréchet derivative of F at x and is also denoted by F′x. Its value at h is denoted by F′xh. Note that F′ is a linear map from X to LX,Y.Definition 2 (see [33, 34]).
F is said to be twice differentiable if the map F′:X⟶LX,Y is Fréchet differentiable. The second derivative of F is denoted by F″ and is a linear map from X to LX,LX,Y that is isomorphic to LX×X,Y.Theorem 3 (see [33, 34]).
The mapF″x∈LX2,Y is symmetric, thatF″xx1,x2=F″xx2,x1, where x1,x2∈X.In this manner,F3x,F4x,⋯ is inductively defined, and Fnx∈LXn,Y is multilinear and symmetric map.Theorem 4 (see [33, 34]).
[Taylor’s theorem]). Suppose thatF∈CnU, where U is an open subset of X containing the line segment from x0 to h, and then(17)Fx0+h=Fx0+F′x0h+12!F″x0h,h+⋯+1n−1!Fn−1x0h,⋯h⏟n−1times+1n−1!∫011−tn−1Fnx0+thh,⋯h⏟ntimesdt+qx,=∑k=0n1n−1!Fn−1x0h,⋯h⏟ktimes,where qx is such that qx=oxn.SinceFkx is symmetric, we denote h,⋯h⏟ktimesbyhk.Using Theorem4,(18)G1=Nu0+u1−Nu0=Nu0+N′u0u1+N″u0u122!+⋯−Nu0=∑k=1∞Nku0u1kk!,(19)G2=Nu0+u1+u2−Nu0+u1=N′u0+u1u2+N″u0+u1u222!+⋯=∑j=1∞∑i=0∞Ni+ju0u1ii!u2jj!,(20)G3=∑i3=1∞∑i2=0∞∑i1=0∞Ni1+i2+i3u0u3i3i3!u2i2i2!u1i1i1!.In general,(21)Gn=∑in=1∞∑in−1=0∞⋯∑i1=0∞N∑k=1niku0∏j=1nujijij!.In view of (9)-(15) and (18)-(21),(22)Nu=G0+G1+G2+G3+⋯=Nu0+∑k=1∞Nku0u1kk!+∑j=1∞∑i=0∞Ni+ju0u1ii!u2jj!+∑i3=1∞∑i2=0∞∑i1=0∞Ni1+i2+i3u0u3i3i3!u2i2i2!u1i1i1!+⋯=Nu0+N′u0u1+u2+u3+⋯+N″u02!u1+u2+u3+⋯2+N3u03!u1+u2+u3+⋯3+⋯.Equation (22) is Taylor series expansion of Nu around u0. Thus, NIM is equivalent to Taylor series expansion around u0.


### 3.3. Error Bound

Thenth degree NIM polynomial at u=u0 is(23)Nnu=Nu0+N′u0u1+u2+u3+⋯+N″u02!u1+u2+u3+⋯2+N3u03!u1+u2+u3+⋯3+⋯+Nnu03!u1+u2+u3+⋯n.Since the NIM approximation becomes more accurate as more terms are included, theNn+1u polynomial must be more accurate than Nnu.(24)Nn+1u=Nu0+N′u01!u1+u2+u3+⋯+N″u02!u1+u2+u3+⋯2+⋯+Nnu0n!u1+u2+u3+⋯n+Nn+1u0n+1!u1+u2+u3+⋯n+1=Nnu+Nn+1u0n+1!u1+u2+u3+⋯n+1.Since the difference betweenNnu and Nn+1u is just that last term, the error of Nnu can be no larger than that term. In other words, the error En is(25)Enu=maxNn+1u0n+1!u1+u2+u3+⋯n+1.Sinceu0 and n are constant in this formula, terms depending only on those constants and u are unaffected by the max operator and can be pulled outside:(26)Enu=maxNn+1u0n+1!u1+u2+u3+⋯n+1.The largest value obtainable byNn+1 could not possibly exceed the maximum value of that derivative between u0 and u. Call the u value that provides that maximum value ω, and the error becomes(27)Enu=Nn+1ωn+1!u1+u2+u3+⋯n+1forallω∈u0,u.


### 3.4. Convergence Analysis

Theorem 5.
Letχ,γ and u be the elements in a Banach Space B, and N is nonlinear contraction from B⟶B such that u=Nχ−Nγ≤kχ−γ, and 0<k<1. Then, according to the principle of Banach fixed point theorem [35], we can prove that(28)um+1≤km+1u0,m=0,1,2,⋯In view of (11)-(14), we have(29)u0=f,u1=Nu0≤ku0⇒u1u0≤k,u2=Nu0+u1−Nu0≤ku1≤k2u0⇒u2u1≤k,u3=Nu0+u1+u2−Nu0+u1≤ku2≤k3u0⇒u3u2≤k,⋮um+1=N∑j=0muj−N∑j=0m−1uj≤kum≤km+1u0⇒um+1um≤k,m=0,1,2,..According to the Banach contraction principle, these are the necessary and sufficient conditions for convergence of the solutions of the NIM procedure. Therefore, the series∑i=0∞um is absolutely and uniformly convergent to the solution of equation (7).


## 3.1. New Iterative Method (NIM)

In this section, we have discussed the new iterative method (NIM) [31, 32] as follows. Let us consider the general equation:(7)uy=Nuy+fy,where N is a nonlinear operator from a Banach space B⟶B, and f is a known function.Lety=y1,y2,⋯..,yn be a solution ofu for equation (7) having the series form.(8)uy=∑i=0∞uiy.The nonlinear operatorN can be written as(9)N∑i=0∞ui=Nu0+∑i=1∞N∑j=0iuj−N∑j=0i−1uj.Then, series solution becomes(10)∑i=1∞ui=f+Nu0+∑i=1∞N∑j=0iuj−N∑j=0i−1uj.We define the consecutive relation:(11)u0=f,(12)u1=Nu0,(13)u2=Nu0+u1−Nu0,(14)um+1=Nu0+⋯+um−Nu0+⋯+um−1,m=1,2,⋯Therefore,(15)u=f+∑i=1∞ui.


## 3.2. Stability Analysis

LetX,Y be Banach spaces and F:X⟶Y a map. LX,Y denotes the set of all linear maps from X to Y. LX,Y is also a Banach space.Definition 1 (see [33, 34]).
F is said to be Fréchet differentiable at x∈X if there exists a continuous linear map A:X⟶Y such that(16)Fx+h−Fx=Ah+wx,h,where limh⟶0wx,h/h=0.A is called the Fréchet derivative of F at x and is also denoted by F′x. Its value at h is denoted by F′xh. Note that F′ is a linear map from X to LX,Y.Definition 2 (see [33, 34]).
F is said to be twice differentiable if the map F′:X⟶LX,Y is Fréchet differentiable. The second derivative of F is denoted by F″ and is a linear map from X to LX,LX,Y that is isomorphic to LX×X,Y.Theorem 3 (see [33, 34]).
The mapF″x∈LX2,Y is symmetric, thatF″xx1,x2=F″xx2,x1, where x1,x2∈X.In this manner,F3x,F4x,⋯ is inductively defined, and Fnx∈LXn,Y is multilinear and symmetric map.Theorem 4 (see [33, 34]).
[Taylor’s theorem]). Suppose thatF∈CnU, where U is an open subset of X containing the line segment from x0 to h, and then(17)Fx0+h=Fx0+F′x0h+12!F″x0h,h+⋯+1n−1!Fn−1x0h,⋯h⏟n−1times+1n−1!∫011−tn−1Fnx0+thh,⋯h⏟ntimesdt+qx,=∑k=0n1n−1!Fn−1x0h,⋯h⏟ktimes,where qx is such that qx=oxn.SinceFkx is symmetric, we denote h,⋯h⏟ktimesbyhk.Using Theorem4,(18)G1=Nu0+u1−Nu0=Nu0+N′u0u1+N″u0u122!+⋯−Nu0=∑k=1∞Nku0u1kk!,(19)G2=Nu0+u1+u2−Nu0+u1=N′u0+u1u2+N″u0+u1u222!+⋯=∑j=1∞∑i=0∞Ni+ju0u1ii!u2jj!,(20)G3=∑i3=1∞∑i2=0∞∑i1=0∞Ni1+i2+i3u0u3i3i3!u2i2i2!u1i1i1!.In general,(21)Gn=∑in=1∞∑in−1=0∞⋯∑i1=0∞N∑k=1niku0∏j=1nujijij!.In view of (9)-(15) and (18)-(21),(22)Nu=G0+G1+G2+G3+⋯=Nu0+∑k=1∞Nku0u1kk!+∑j=1∞∑i=0∞Ni+ju0u1ii!u2jj!+∑i3=1∞∑i2=0∞∑i1=0∞Ni1+i2+i3u0u3i3i3!u2i2i2!u1i1i1!+⋯=Nu0+N′u0u1+u2+u3+⋯+N″u02!u1+u2+u3+⋯2+N3u03!u1+u2+u3+⋯3+⋯.Equation (22) is Taylor series expansion of Nu around u0. Thus, NIM is equivalent to Taylor series expansion around u0.


## 3.3. Error Bound

Thenth degree NIM polynomial at u=u0 is(23)Nnu=Nu0+N′u0u1+u2+u3+⋯+N″u02!u1+u2+u3+⋯2+N3u03!u1+u2+u3+⋯3+⋯+Nnu03!u1+u2+u3+⋯n.Since the NIM approximation becomes more accurate as more terms are included, theNn+1u polynomial must be more accurate than Nnu.(24)Nn+1u=Nu0+N′u01!u1+u2+u3+⋯+N″u02!u1+u2+u3+⋯2+⋯+Nnu0n!u1+u2+u3+⋯n+Nn+1u0n+1!u1+u2+u3+⋯n+1=Nnu+Nn+1u0n+1!u1+u2+u3+⋯n+1.Since the difference betweenNnu and Nn+1u is just that last term, the error of Nnu can be no larger than that term. In other words, the error En is(25)Enu=maxNn+1u0n+1!u1+u2+u3+⋯n+1.Sinceu0 and n are constant in this formula, terms depending only on those constants and u are unaffected by the max operator and can be pulled outside:(26)Enu=maxNn+1u0n+1!u1+u2+u3+⋯n+1.The largest value obtainable byNn+1 could not possibly exceed the maximum value of that derivative between u0 and u. Call the u value that provides that maximum value ω, and the error becomes(27)Enu=Nn+1ωn+1!u1+u2+u3+⋯n+1forallω∈u0,u.


## 3.4. Convergence Analysis

Theorem 5.
Letχ,γ and u be the elements in a Banach Space B, and N is nonlinear contraction from B⟶B such that u=Nχ−Nγ≤kχ−γ, and 0<k<1. Then, according to the principle of Banach fixed point theorem [35], we can prove that(28)um+1≤km+1u0,m=0,1,2,⋯In view of (11)-(14), we have(29)u0=f,u1=Nu0≤ku0⇒u1u0≤k,u2=Nu0+u1−Nu0≤ku1≤k2u0⇒u2u1≤k,u3=Nu0+u1+u2−Nu0+u1≤ku2≤k3u0⇒u3u2≤k,⋮um+1=N∑j=0muj−N∑j=0m−1uj≤kum≤km+1u0⇒um+1um≤k,m=0,1,2,..According to the Banach contraction principle, these are the necessary and sufficient conditions for convergence of the solutions of the NIM procedure. Therefore, the series∑i=0∞um is absolutely and uniformly convergent to the solution of equation (7).


## 4. Application of NIM

Here, we monitor the pollutants for three different types of input models. According to the values ofqt, these three models are impulse input, step input, and sinusoidal input [14, 30];


### 4.1. Impulse Input

This model is considered when there is a spike in concentration of dumped contaminant. Considering the value ofqt is 100 units, the system becomes(30)dz1dt=381180z3t+100−382900z1t,(31)dz2dt=182900z1t−18850z2t,(32)dz3dt=202900z1t+18850z2t−381180z3t,(33)withinitialconditionsz10=z20=z30=0.Integrating the system of equations from 0 tot and using the initial conditions (33), we get(34)z1t=100t+∫0t19z3t590dt−∫0t19z1t1450dt,z2t=∫0t9z1t1450dt−∫0t9z2t425dt,z3t=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.We let the nonlinear terms as(35)N1z1,z2,z3=∫0t19z3t590dt−∫0t19z1t1450dt,N2z1,z2,z3=∫0t9z1t1450dt−∫0t9z2t425dt,N3z1,z2,z3=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.In view of (11)–(14), we obtain first few components of the new iterative solution for the considered problem as follows:(36)z1,0t=100t,z2,0t=0,z3,0t=0,z1,1t=N1z1,0,z2,0,z3,0=∫0t19z3,0t590dt−∫0t19z1,0t1450dt=−19t229,z2,1t=N2z1,0,z2,0,z3,0=∫0t9z1,0t1450dt−∫0t9z2,0t425dt=9t229,z3,1t=N3z1,0,z2,0,z3,0=∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=10t229,z1,2t=N1z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t19z3,0t+z3,1t590dt−∫0t19z1,0t+z1,1t1450dt−∫0t19z3,0t590dt−∫0t19z1,0t1450dt=16283t32480950,z2,2t=N2z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t9z1,0t+z1,1t1450dt−∫0t9z2,0t+z2,1t425dt−∫0t9z1,0t1450dt−∫0t9z2,0t425dt=−507t3142970,z3,2t=N3z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0tz1,0t+z1,1t145dt+∫0t9z2,0t+z2,1t425dt−∫0t19z3,0t+z3,1t590dt−∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=−63623t321088075and so on. In the same manner, the successive approximations can be obtained.Therefore, the 9-iteration approximate series solutions are(37)z1t=∑i=09z1,it=100t−0.6551724138t2+0.006563211673t3−0.00004578968898t4+2.284062102×10−7t5−7.613669180×10−10t6+9.264952117×10−13t7+8.062395960×10−15t8−7.447533873×10−17t9+3.912821015×10−19t10,z2t=∑i=09z2,it=0.3103448276t2−0.003546198503t3+0.00002895828605t4−1.794892311×10−7t5+8.697736903×10−10t6−3.306351810×10−12t7+9.470940223×10−15t8−1.672429216×10−17t9−1.080992422×10−20t10,z3t=∑i=09z3,i=0.3448275862t2−0.003017013170t3+0.00001683140293t4−4.891697909×10−8t5−1.084067722×10−10t6+2.379856598×10−12t7−1.753333618×10−14t8+9.119963089×10−17t9−3.804721773×10−19t10.It is easily observed that the obtained solutionz1,z2,z3=∑i=0∞z1,i,∑i=0∞z2,i,∑i=0∞z3,i, converges to the exact solution when(38)β1,i=z1,i+1z1,i<1,β2,i=z2,i+1z2,i<1,β3,i=z3,i+1z3,i<1,fori≥1,0<t≤1.Here, we can see that(39)z1,2z1,1=16283t3/248095019t2/29=0.0100175336<1,z2,2z2,1=507t3/1429709t2/29=0.01142663961<1,z3,2z3,1=63623t3/2108807510t2/29=0.00874933816<1.This confirm that the NIM series solution converges to the exact solution.


### 4.2. Step Input

This model is considered when contaminants enter the lake continuously at steady concentration. Considering the value ofqt is 100 t units, the system becomes(40)dz1dt=381180z3t+100t−382900z1t,dz2dt=182900z1t−18850z2t,dz3dt=202900z1t+18850z2t−381180z3t,with initial conditions(41)z10=z20=z30=0.Integrating the system of equations from 0 tot and using the initial conditions (41), we get(42)z1t=50t2+∫0t19z3t590dt−∫0t19z1t1450dt,z2t=∫0t9z1t1450dt−∫0t9z2t425dt,z3t=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.We let the nonlinear terms as(43)N1z1,z2,z3=∫0t19z3t590dt−∫0t19z1t1450dt,N2z1,z2,z3=∫0t9z1t1450dt−∫0t9z2t425dt,N3z1,z2,z3=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.In view of (11)–(14), we have the following first few components of the new iterative solution for the considered problem:(44)z1,0t=50t2,z2,0t=0,z3,0t=0,z1,1t=N1z1,0,z2,0,z3,0=∫0t19z3,0t590dt−∫0t19z1,0t1450dt=−19t287,z2,1t=N2z1,0,z2,0,z3,0=∫0t9z1,0t1450dt−∫0t9z2,0t425dt=3t229,z3,1t=N3z1,0,z2,0,z3,0=∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=10t287,z1,2t=N1z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t19z3,0t+z3,1t590dt−∫0t19z1,0t+z1,1t1450dt−∫0t19z3,0t590dt−∫0t19z1,0t1450dt=16283t39923800,z2,2t=N2z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t9z1,0t+z1,1t1450dt−∫0t9z2,0t+z2,1t425dt−∫0t9z1,0t1450dt−∫0t9z2,0t425dt=−507t3571880,z3,2t=N3z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0tz1,0t+z1,1t145dt+∫0t9z2,0t+z2,1t425dt−∫0t19z3,0t+z3,1t590dt−∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=−63623t384352300,and so on. In the same manner, the successive components can be obtained.Therefore, the 9-iteration approximate series solutions are(45)z1t=∑i=09z1,it=50t2−0.2183908046t3+0.001640802918t4−9.157937795×10−6t5+3.80677017×10−10t6−1.087667026×10−10t7+1.158119015×10−13t8+8.958217733×10−16t9−7.447533873×10−18t10+3.557110013×10−20t11,z2t=∑i=09z2,it=0.1034482759t3−0.0008865496258t4+5.791657210×10−6t5−2.991487185×10−8t6+1.242533843×10−10t7−4.132939762×10−13t8+1.052326691×10−15t9−1.672429216×10−18t10−9.827203838×10−22t11,z3t=∑i=09z3,i=0.1149425287t3−0.0007542532924t4+3.366280585×10−6t5–8.152829849×10−9t6–1.548668175×10−11t7+2.974820747×10−13t8–1.948148465×10−15t9+9.119963089×10−18t10–3.458837975×10−20t11.It is easily observed that the obtained solution,z1,z2,z3=∑i=0∞z1,i,∑i=0∞z2,i,∑i=0∞z3,i, converges to the exact solution when(46)β1,i=z1,i+1z1,i<1,β2,i=z2,i+1z2,i<1,β3,i=z3,i+1z3,i<1,fori≥1,0<t≤1.Here, we can see that(47)z1,2z1,1=16283t3/9923800−19t2/87=0.00751315016<1,z2,2z2,1=−507t3/5718803t2/29=0.02916446054<1,z3,2z3,1=−63623t3/8435230010t2/87=0.00656200362<1.This confirms that the NIM series solution converges to the exact solution.


### 4.3. Sinusoidal Input

When concentration of the lake converges to the average input concentration of the contaminant [14], then the input is periodic input, and we assume qt=1+sint(48)dz1dt=381180z3t+1+sint−382900z1t,(49)dz2dt=182900z1t−18850z2t,(50)dz3dt=202900z1t+18850z2t−381180z3t,(51)withinitialconditionsz10=z20=z30=0.Integrating the system of equations from 0 tot and using the initial conditions (51), we get(52)z1t=1+t−cost+∫0t19z3t590dt−∫0t19z1t1450dt,z2t=∫0t9z1t1450dt−∫0t9z2t425dt,z3t=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.We let the nonlinear terms as(53)N1z1,z2,z3=∫0t19z3t590dt−∫0t19z1t1450dt,N2z1,z2,z3=∫0t9z1t1450dt−∫0t9z2t425dt,N3z1,z2,z3=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.In view of (11)–(14), we have the following first few components of the new iterative solution for the considered problem:(54)z1,0t=1+t−cost,z2,0t=0,(63)(55)z3,0t=0,z1,1t=N1z1,0,z2,0,z3,0=∫0t19z3,0t590dt−∫0t19z1,0t1450dt=−19t1450−19t22900+19sint1450,z2,1t=N2z1,0,z2,0,z3,0=∫0t9z1,0t1450dt−∫0t9z2,0t425dt=9t1450+9t22900−9sint1450,z3,1t=N3z1,0,z2,0,z3,0=∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=t145+t2290−sint145,z1,2t=N1z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t19z3,0t+z3,1t590dt−∫0t19z1,0t+z1,1t1450dt−∫0t19z3,0t590dt−∫0t19z1,0t1450dt=−48849124047500+48849t2248095000+16283t3248095000+48849cost124047500,z2,2t=N2z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t9z1,0t+z1,1t1450dt−∫0t9z2,0t+z2,1t425dt−∫0t9z1,0t1450dt−∫0t9z2,0t425dt=15217148500−1521t214297000−507t314297000−1521cost7148500,z3,2t=N3z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0tz1,0t+z1,1t145dt+∫0t9z2,0t+z2,1t425dt−∫0t19z3,0t+z3,1t590dt−∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=1908691054403750−190869t22108807500−63623t32108807500−190869cost1054403750,and so on. In the same manner, the rest of components can be obtained.Therefore, the 9-iteration approximate series solutions are(56)z1t=∑i=09z1,it=0.01309246424sint−7.548439160×10−12t6+0.9869075358t+7.392117911×10−17t8+2.238353115×10−9t5−4.464785362×10−7t4+9.907126563×10−15t7+0.00006380144335×t3−0.006354964808×t2−0.9996064813cost+0.9996064813+3.912821015×10−21t10−7.056251772×10−19t9,z2t=∑i=09z2,it=−0.006199952822sint+8.466376566×10−12t6+0.006199952822t+9.320421594×10−17t8–1.742737715×10−9t5+2.806153398×10−7t4−3.230576505×10−14t7+0.0002125566897–0.00003430469668t3+0.002997169931t2−0.0002125566897cost–1.080992422×10−22t10–1.683239140×10−19t9,z3t=∑i=09z3,i=−0.006892511414sint–9.179374067×10−13t6+0.006892511414t–1.671253951×10−16t8−4.956154003×10−10t5+1.658631964×10−7t4+2.239863848×10−14t7−0.00002949674667t3+0.003357794877t2–0.0001809619702cost–3.804721773×10−21t10+8.739490912×10−19t9+0.0001809619702.It is easily observed that the obtained solution,z1,z2,z3=∑i=0∞z1,i,∑i=0∞z2,i,∑i=0∞z3,i, converges to the exact solution when(57)β1,i=z1,i+1z1,i<1,β2,i=z2,i+1z2,i<1,β3,i=z3,i+1z3,i<1,fori≥1,0<t≤1.Here, we can see that(58)z1,2z1,1=−48849/124047500+48849t2/248095000+16283t3/248095000+48849cost/124047500−19t/1450−19t2/2900+19sint/1450=0.00944522675<1,z2,2z2,1=1521/7148500−1521t2/14297000−507t3/14297000−1521cost/71485009t/1450+9t2/2900−9sint/1450=0.010773829<1,z3,2z3,1=190869/1054403750−190869t2/2108807500−63623t3/2108807500−190869cost/1054403750t/145+t2/290−sint/145=0.00824948252<1.This confirms that the NIM series solution converges to the exact solution.


## 4.1. Impulse Input

This model is considered when there is a spike in concentration of dumped contaminant. Considering the value ofqt is 100 units, the system becomes(30)dz1dt=381180z3t+100−382900z1t,(31)dz2dt=182900z1t−18850z2t,(32)dz3dt=202900z1t+18850z2t−381180z3t,(33)withinitialconditionsz10=z20=z30=0.Integrating the system of equations from 0 tot and using the initial conditions (33), we get(34)z1t=100t+∫0t19z3t590dt−∫0t19z1t1450dt,z2t=∫0t9z1t1450dt−∫0t9z2t425dt,z3t=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.We let the nonlinear terms as(35)N1z1,z2,z3=∫0t19z3t590dt−∫0t19z1t1450dt,N2z1,z2,z3=∫0t9z1t1450dt−∫0t9z2t425dt,N3z1,z2,z3=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.In view of (11)–(14), we obtain first few components of the new iterative solution for the considered problem as follows:(36)z1,0t=100t,z2,0t=0,z3,0t=0,z1,1t=N1z1,0,z2,0,z3,0=∫0t19z3,0t590dt−∫0t19z1,0t1450dt=−19t229,z2,1t=N2z1,0,z2,0,z3,0=∫0t9z1,0t1450dt−∫0t9z2,0t425dt=9t229,z3,1t=N3z1,0,z2,0,z3,0=∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=10t229,z1,2t=N1z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t19z3,0t+z3,1t590dt−∫0t19z1,0t+z1,1t1450dt−∫0t19z3,0t590dt−∫0t19z1,0t1450dt=16283t32480950,z2,2t=N2z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t9z1,0t+z1,1t1450dt−∫0t9z2,0t+z2,1t425dt−∫0t9z1,0t1450dt−∫0t9z2,0t425dt=−507t3142970,z3,2t=N3z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0tz1,0t+z1,1t145dt+∫0t9z2,0t+z2,1t425dt−∫0t19z3,0t+z3,1t590dt−∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=−63623t321088075and so on. In the same manner, the successive approximations can be obtained.Therefore, the 9-iteration approximate series solutions are(37)z1t=∑i=09z1,it=100t−0.6551724138t2+0.006563211673t3−0.00004578968898t4+2.284062102×10−7t5−7.613669180×10−10t6+9.264952117×10−13t7+8.062395960×10−15t8−7.447533873×10−17t9+3.912821015×10−19t10,z2t=∑i=09z2,it=0.3103448276t2−0.003546198503t3+0.00002895828605t4−1.794892311×10−7t5+8.697736903×10−10t6−3.306351810×10−12t7+9.470940223×10−15t8−1.672429216×10−17t9−1.080992422×10−20t10,z3t=∑i=09z3,i=0.3448275862t2−0.003017013170t3+0.00001683140293t4−4.891697909×10−8t5−1.084067722×10−10t6+2.379856598×10−12t7−1.753333618×10−14t8+9.119963089×10−17t9−3.804721773×10−19t10.It is easily observed that the obtained solutionz1,z2,z3=∑i=0∞z1,i,∑i=0∞z2,i,∑i=0∞z3,i, converges to the exact solution when(38)β1,i=z1,i+1z1,i<1,β2,i=z2,i+1z2,i<1,β3,i=z3,i+1z3,i<1,fori≥1,0<t≤1.Here, we can see that(39)z1,2z1,1=16283t3/248095019t2/29=0.0100175336<1,z2,2z2,1=507t3/1429709t2/29=0.01142663961<1,z3,2z3,1=63623t3/2108807510t2/29=0.00874933816<1.This confirm that the NIM series solution converges to the exact solution.


## 4.2. Step Input

This model is considered when contaminants enter the lake continuously at steady concentration. Considering the value ofqt is 100 t units, the system becomes(40)dz1dt=381180z3t+100t−382900z1t,dz2dt=182900z1t−18850z2t,dz3dt=202900z1t+18850z2t−381180z3t,with initial conditions(41)z10=z20=z30=0.Integrating the system of equations from 0 tot and using the initial conditions (41), we get(42)z1t=50t2+∫0t19z3t590dt−∫0t19z1t1450dt,z2t=∫0t9z1t1450dt−∫0t9z2t425dt,z3t=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.We let the nonlinear terms as(43)N1z1,z2,z3=∫0t19z3t590dt−∫0t19z1t1450dt,N2z1,z2,z3=∫0t9z1t1450dt−∫0t9z2t425dt,N3z1,z2,z3=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.In view of (11)–(14), we have the following first few components of the new iterative solution for the considered problem:(44)z1,0t=50t2,z2,0t=0,z3,0t=0,z1,1t=N1z1,0,z2,0,z3,0=∫0t19z3,0t590dt−∫0t19z1,0t1450dt=−19t287,z2,1t=N2z1,0,z2,0,z3,0=∫0t9z1,0t1450dt−∫0t9z2,0t425dt=3t229,z3,1t=N3z1,0,z2,0,z3,0=∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=10t287,z1,2t=N1z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t19z3,0t+z3,1t590dt−∫0t19z1,0t+z1,1t1450dt−∫0t19z3,0t590dt−∫0t19z1,0t1450dt=16283t39923800,z2,2t=N2z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t9z1,0t+z1,1t1450dt−∫0t9z2,0t+z2,1t425dt−∫0t9z1,0t1450dt−∫0t9z2,0t425dt=−507t3571880,z3,2t=N3z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0tz1,0t+z1,1t145dt+∫0t9z2,0t+z2,1t425dt−∫0t19z3,0t+z3,1t590dt−∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=−63623t384352300,and so on. In the same manner, the successive components can be obtained.Therefore, the 9-iteration approximate series solutions are(45)z1t=∑i=09z1,it=50t2−0.2183908046t3+0.001640802918t4−9.157937795×10−6t5+3.80677017×10−10t6−1.087667026×10−10t7+1.158119015×10−13t8+8.958217733×10−16t9−7.447533873×10−18t10+3.557110013×10−20t11,z2t=∑i=09z2,it=0.1034482759t3−0.0008865496258t4+5.791657210×10−6t5−2.991487185×10−8t6+1.242533843×10−10t7−4.132939762×10−13t8+1.052326691×10−15t9−1.672429216×10−18t10−9.827203838×10−22t11,z3t=∑i=09z3,i=0.1149425287t3−0.0007542532924t4+3.366280585×10−6t5–8.152829849×10−9t6–1.548668175×10−11t7+2.974820747×10−13t8–1.948148465×10−15t9+9.119963089×10−18t10–3.458837975×10−20t11.It is easily observed that the obtained solution,z1,z2,z3=∑i=0∞z1,i,∑i=0∞z2,i,∑i=0∞z3,i, converges to the exact solution when(46)β1,i=z1,i+1z1,i<1,β2,i=z2,i+1z2,i<1,β3,i=z3,i+1z3,i<1,fori≥1,0<t≤1.Here, we can see that(47)z1,2z1,1=16283t3/9923800−19t2/87=0.00751315016<1,z2,2z2,1=−507t3/5718803t2/29=0.02916446054<1,z3,2z3,1=−63623t3/8435230010t2/87=0.00656200362<1.This confirms that the NIM series solution converges to the exact solution.


## 4.3. Sinusoidal Input

When concentration of the lake converges to the average input concentration of the contaminant [14], then the input is periodic input, and we assume qt=1+sint(48)dz1dt=381180z3t+1+sint−382900z1t,(49)dz2dt=182900z1t−18850z2t,(50)dz3dt=202900z1t+18850z2t−381180z3t,(51)withinitialconditionsz10=z20=z30=0.Integrating the system of equations from 0 tot and using the initial conditions (51), we get(52)z1t=1+t−cost+∫0t19z3t590dt−∫0t19z1t1450dt,z2t=∫0t9z1t1450dt−∫0t9z2t425dt,z3t=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.We let the nonlinear terms as(53)N1z1,z2,z3=∫0t19z3t590dt−∫0t19z1t1450dt,N2z1,z2,z3=∫0t9z1t1450dt−∫0t9z2t425dt,N3z1,z2,z3=∫0tz1t145dt+∫0t9z2t425dt−∫0t19z3t590dt.In view of (11)–(14), we have the following first few components of the new iterative solution for the considered problem:(54)z1,0t=1+t−cost,z2,0t=0,(63)(55)z3,0t=0,z1,1t=N1z1,0,z2,0,z3,0=∫0t19z3,0t590dt−∫0t19z1,0t1450dt=−19t1450−19t22900+19sint1450,z2,1t=N2z1,0,z2,0,z3,0=∫0t9z1,0t1450dt−∫0t9z2,0t425dt=9t1450+9t22900−9sint1450,z3,1t=N3z1,0,z2,0,z3,0=∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=t145+t2290−sint145,z1,2t=N1z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t19z3,0t+z3,1t590dt−∫0t19z1,0t+z1,1t1450dt−∫0t19z3,0t590dt−∫0t19z1,0t1450dt=−48849124047500+48849t2248095000+16283t3248095000+48849cost124047500,z2,2t=N2z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0t9z1,0t+z1,1t1450dt−∫0t9z2,0t+z2,1t425dt−∫0t9z1,0t1450dt−∫0t9z2,0t425dt=15217148500−1521t214297000−507t314297000−1521cost7148500,z3,2t=N3z1,0+z1,1,z2,0+z2,1,z3,0+z3,1−N1z1,0,z2,0,z3,0=∫0tz1,0t+z1,1t145dt+∫0t9z2,0t+z2,1t425dt−∫0t19z3,0t+z3,1t590dt−∫0tz1,0t145dt+∫0t9z2,0t425dt−∫0t19z3,0t590dt=1908691054403750−190869t22108807500−63623t32108807500−190869cost1054403750,and so on. In the same manner, the rest of components can be obtained.Therefore, the 9-iteration approximate series solutions are(56)z1t=∑i=09z1,it=0.01309246424sint−7.548439160×10−12t6+0.9869075358t+7.392117911×10−17t8+2.238353115×10−9t5−4.464785362×10−7t4+9.907126563×10−15t7+0.00006380144335×t3−0.006354964808×t2−0.9996064813cost+0.9996064813+3.912821015×10−21t10−7.056251772×10−19t9,z2t=∑i=09z2,it=−0.006199952822sint+8.466376566×10−12t6+0.006199952822t+9.320421594×10−17t8–1.742737715×10−9t5+2.806153398×10−7t4−3.230576505×10−14t7+0.0002125566897–0.00003430469668t3+0.002997169931t2−0.0002125566897cost–1.080992422×10−22t10–1.683239140×10−19t9,z3t=∑i=09z3,i=−0.006892511414sint–9.179374067×10−13t6+0.006892511414t–1.671253951×10−16t8−4.956154003×10−10t5+1.658631964×10−7t4+2.239863848×10−14t7−0.00002949674667t3+0.003357794877t2–0.0001809619702cost–3.804721773×10−21t10+8.739490912×10−19t9+0.0001809619702.It is easily observed that the obtained solution,z1,z2,z3=∑i=0∞z1,i,∑i=0∞z2,i,∑i=0∞z3,i, converges to the exact solution when(57)β1,i=z1,i+1z1,i<1,β2,i=z2,i+1z2,i<1,β3,i=z3,i+1z3,i<1,fori≥1,0<t≤1.Here, we can see that(58)z1,2z1,1=−48849/124047500+48849t2/248095000+16283t3/248095000+48849cost/124047500−19t/1450−19t2/2900+19sint/1450=0.00944522675<1,z2,2z2,1=1521/7148500−1521t2/14297000−507t3/14297000−1521cost/71485009t/1450+9t2/2900−9sint/1450=0.010773829<1,z3,2z3,1=190869/1054403750−190869t2/2108807500−63623t3/2108807500−190869cost/1054403750t/145+t2/290−sint/145=0.00824948252<1.This confirms that the NIM series solution converges to the exact solution.


## 5. Results and Discussions

In this work, we have presented the new iterative method (NIM) as a useful analytic-numeric tool to solve pollution model for a system of three interconnecting lakes. The NIM algorithm is coded in the computer algebra package Maple software. The Maple condition variable digits controlling the number of significant digits are set to 16 in all the calculations done in this paper. We additionally fix the estimations of the parameters and take the initial conditions as stated in [14]. The time domain studied in this work is [0, 20], same as considered in [20]. For comparison purposes, the fourth-order Runge-Kutta algorithm (RK4) build in Maple software was used to obtain the numeric solution of the pollution problems.


### 5.1. Impulse Input

In Table1, we have made tabular comparison between numeric solutions of NIM and ADM for the impulse input model for t up to 0.1 only to directly compare with [17]. The 9-iteration NIM solutions agree with the RK4 solutions at least 7 decimal places as shown in Table 2 for impulse input for t up to 20. We have figured out the z1-z2, z2-z3, z1-z3, and z1-z2-z3 phase portraits of 9-iteration NIM for the impulse input system in Figure 2 which shows their smooth behavior in the phase plane. Figure 3 shows the graphical comparison among RK4 solutions and the NIM approximations for the impulse input models for t up to 20.Table 1
Numerical comparisons between ADM and NIM for the impulse input model.
ADM [14]NIMtz1tz2tz3tz1tz2tz3t0.010.999930.000030.000030.999930.000030.000030.021.999740.000120.000141.999740.000120.000140.032.999410.000280.000312.999410.000280.000310.043.999000.000490.000553.998950.000500.000550.054.998950.000760.000864.998360.000780.000860.065.998360.001120.001245.997640.001120.001240.076.997640.001520.001696.996790.001520.001690.087.996790.001980.002217.995810.001980.002210.098.995810.002510.002798.994690.002510.002790.19.993450.003100.003459.993450.003100.00345Table 2
Differences between 9-iteration NIM with RK4 solutions in the impulse input model.
∆=RK40.001−NIMt∆z1t∆z2t∆z3t21.02e-111.28e-131.42e-1343.02e-117.62e-138.56e-1363.19e-111.08e-128.2e-1381.038e-102.39e-126.27e-12101.435e-102.349e-111.306e-10121.028e-091.763e-109.677e-10145.914e-099.273e-105.216e-09162.608e-083.946e-092.248e-08189.562e-081.419e-088.153e-08203.031e-074.472e-082.582e-07Figure 2
Phase portraits using 9-iteration NIM for impulse input.Figure 3
Comparison among NIM (solid line) and RK4 (circle) for (a)z1t, (b) z2t, and (c) z3tin the impulse input model.
(a)(b)(c)


### 5.2. Step Input

In Table3, we have made tabular comparison between numeric solutions of NIM and ADM for the step input model for t up to 0.1 only as [17]. We have noted from Table 4 that the NIM precisions match to those of RK4 at least 6 decimal places for t up to 20. The routine was configured to find out the graphical representation of absolute errors (A.E.) for the step input model as it can be seen in Figure 4. In general terms, the absolute error is low for all approximations. If required, the error can be further reduced by increasing the iterations of the NIM approximations. Figure 5 shows the graphical comparison among RK4 solutions and the NIM approximations for the step input models for t up to 20.Table 3
Numerical comparisons between ADM and NIM for the step input model.
ADM [14]NIMtz1tz2tz3tz1tz2tz3t0.010.005000.10344×10−60.11493×10−60.005000.10344×10−60.11493×10−60.020.020000.82744×10−60.91942×10−60.020000.82744×10−60.91942×10−60.030.044990.27923×10−50.31028×10−50.044990.2792×10−50.31028×10−50.040.079990.66184×10−50.73545×10−50.079990.66184×10−50.73544×10−50.050.124970.000010.000010.124970.000010.000010.060.179950.000020.000020.179950.000020.000020.070.244930.000040.000040.244930.000040.000040.080.319890.000050.000060.319890.000050.000060.090.404840.000080.000080.404840.000080.000080.10.499780.000100.000110.499780.000100.00011Table 4
Differences between 9-iteration NIM with RK4 solutions in the step input model.
∆=∣RK40.001−NIM∣t∆z1t∆z2t∆z3t21.60e-111.203e-131.365e-1348.11e-111.265e-121.421e-1262.45e-104.48e-124.88e-1286.81e-101.780e-111.424e-11101.50e-102.311e-111.127e-10122.81e-102.196e-101.019e-09144.476e-091.191e-096.223e-09163.119e-085.482e-093.029e-08181.439e-072.148e-081.228e-07205.113e-077.496e-084.322e-07Figure 4
Absolute error (A.E.) of step input for∆z1t, ∆z2t, and ∆z3t.Figure 5
Comparison among NIM (solid line) and RK4 (circle) for (a)z1t, (b) z2t, and (c) z3tin the step input model.
(a)(b)(c)


### 5.3. Sinusoidal Input

In Table5, we have made tabular comparison between numeric solutions of NIM and ADM [17] for the sinusoidal input model for t up to 0.1 same as [17]. Similarly, we have noticed from Table 6 that the NIM precisions match to those of RK4 at least minimum 5 decimal places for sinusoidal inputs for t up to 20. We reproduced the z1-z2, z2-z3, z1-z3, and z1-z2-z3 phase portraits of 9-iteration NIM for the sinusoidal input system in Figure 6 which shows their smooth behavior in the phase plane. Figure 7 shows the graphical comparison among RK4 solutions and the NIM approximations for the sinusoidal input models for t up to 20.Table 5
Numerical comparisons between ADM and NIM for the sinusoidal input model.
ADM [14]NIMtz1tz2tz3tz1tz2tz3t0.010.010050.31134×10−60.34595×10−60.010050.31134×10−60.34595×10−60.020.020200.12494×10−50.13883×10−50.020200.12494×10−50.13883×10−50.030.030440.28201×10−50.31337×10−50.030440.28201×10−50.31337×10−50.040.040790.50294×10−50.55888×10−50.040790.50294×10−50.55888×10−50.050.051230.78834×10−50.87605×10−50.051230.78834×10−50.87605×10−50.060.061780.000010.000010.061780.000010.000010.070.072420.000020.000020.072420.000020.000020.080.083160.000020.000020.083160.000020.000020.090.093990.000020.000030.093990.000030.000030.10.104930.000030.000040.104930.000030.00004Table 6
Differences between 9-iteration NIM with RK4 solutions in the sinusoidal input model.
∆=∣RK40.001−NIM∣t∆z1t∆z2t∆z3t21.3e-1402e-1649e-152e-165e-1661.0e-142e-168.4e-1581.85e-132.77e-141.690e-13102.33e-123.479e-131.974e-12121.726e-112.572e-121.462e-11149.251e-111.378e-117.863e-11163.933e-105.847e-113.349e-10181.404e-092.079e-101.196e-09204.366e-096.444e-103.722e-09Figure 6
Phase portraits using 9-iteration NIM for sinusoidal input.Figure 7
Comparison among 9-iteration NIM (solid line) and RK4 (circle) for (a)z1t, (b) z2t, and (c) z3tin the sinusoidal input model.
(a)(b)(c)All of the above discussions indicate that NIM is very accurate, requires less computation work, very easy to implement, and consumes less time to obtain good approximate solution in a wide interval for lake pollution problem.


## 5.1. Impulse Input

In Table1, we have made tabular comparison between numeric solutions of NIM and ADM for the impulse input model for t up to 0.1 only to directly compare with [17]. The 9-iteration NIM solutions agree with the RK4 solutions at least 7 decimal places as shown in Table 2 for impulse input for t up to 20. We have figured out the z1-z2, z2-z3, z1-z3, and z1-z2-z3 phase portraits of 9-iteration NIM for the impulse input system in Figure 2 which shows their smooth behavior in the phase plane. Figure 3 shows the graphical comparison among RK4 solutions and the NIM approximations for the impulse input models for t up to 20.Table 1
Numerical comparisons between ADM and NIM for the impulse input model.
ADM [14]NIMtz1tz2tz3tz1tz2tz3t0.010.999930.000030.000030.999930.000030.000030.021.999740.000120.000141.999740.000120.000140.032.999410.000280.000312.999410.000280.000310.043.999000.000490.000553.998950.000500.000550.054.998950.000760.000864.998360.000780.000860.065.998360.001120.001245.997640.001120.001240.076.997640.001520.001696.996790.001520.001690.087.996790.001980.002217.995810.001980.002210.098.995810.002510.002798.994690.002510.002790.19.993450.003100.003459.993450.003100.00345Table 2
Differences between 9-iteration NIM with RK4 solutions in the impulse input model.
∆=RK40.001−NIMt∆z1t∆z2t∆z3t21.02e-111.28e-131.42e-1343.02e-117.62e-138.56e-1363.19e-111.08e-128.2e-1381.038e-102.39e-126.27e-12101.435e-102.349e-111.306e-10121.028e-091.763e-109.677e-10145.914e-099.273e-105.216e-09162.608e-083.946e-092.248e-08189.562e-081.419e-088.153e-08203.031e-074.472e-082.582e-07Figure 2
Phase portraits using 9-iteration NIM for impulse input.Figure 3
Comparison among NIM (solid line) and RK4 (circle) for (a)z1t, (b) z2t, and (c) z3tin the impulse input model.
(a)(b)(c)


## 5.2. Step Input

In Table3, we have made tabular comparison between numeric solutions of NIM and ADM for the step input model for t up to 0.1 only as [17]. We have noted from Table 4 that the NIM precisions match to those of RK4 at least 6 decimal places for t up to 20. The routine was configured to find out the graphical representation of absolute errors (A.E.) for the step input model as it can be seen in Figure 4. In general terms, the absolute error is low for all approximations. If required, the error can be further reduced by increasing the iterations of the NIM approximations. Figure 5 shows the graphical comparison among RK4 solutions and the NIM approximations for the step input models for t up to 20.Table 3
Numerical comparisons between ADM and NIM for the step input model.
ADM [14]NIMtz1tz2tz3tz1tz2tz3t0.010.005000.10344×10−60.11493×10−60.005000.10344×10−60.11493×10−60.020.020000.82744×10−60.91942×10−60.020000.82744×10−60.91942×10−60.030.044990.27923×10−50.31028×10−50.044990.2792×10−50.31028×10−50.040.079990.66184×10−50.73545×10−50.079990.66184×10−50.73544×10−50.050.124970.000010.000010.124970.000010.000010.060.179950.000020.000020.179950.000020.000020.070.244930.000040.000040.244930.000040.000040.080.319890.000050.000060.319890.000050.000060.090.404840.000080.000080.404840.000080.000080.10.499780.000100.000110.499780.000100.00011Table 4
Differences between 9-iteration NIM with RK4 solutions in the step input model.
∆=∣RK40.001−NIM∣t∆z1t∆z2t∆z3t21.60e-111.203e-131.365e-1348.11e-111.265e-121.421e-1262.45e-104.48e-124.88e-1286.81e-101.780e-111.424e-11101.50e-102.311e-111.127e-10122.81e-102.196e-101.019e-09144.476e-091.191e-096.223e-09163.119e-085.482e-093.029e-08181.439e-072.148e-081.228e-07205.113e-077.496e-084.322e-07Figure 4
Absolute error (A.E.) of step input for∆z1t, ∆z2t, and ∆z3t.Figure 5
Comparison among NIM (solid line) and RK4 (circle) for (a)z1t, (b) z2t, and (c) z3tin the step input model.
(a)(b)(c)


## 5.3. Sinusoidal Input

In Table5, we have made tabular comparison between numeric solutions of NIM and ADM [17] for the sinusoidal input model for t up to 0.1 same as [17]. Similarly, we have noticed from Table 6 that the NIM precisions match to those of RK4 at least minimum 5 decimal places for sinusoidal inputs for t up to 20. We reproduced the z1-z2, z2-z3, z1-z3, and z1-z2-z3 phase portraits of 9-iteration NIM for the sinusoidal input system in Figure 6 which shows their smooth behavior in the phase plane. Figure 7 shows the graphical comparison among RK4 solutions and the NIM approximations for the sinusoidal input models for t up to 20.Table 5
Numerical comparisons between ADM and NIM for the sinusoidal input model.
ADM [14]NIMtz1tz2tz3tz1tz2tz3t0.010.010050.31134×10−60.34595×10−60.010050.31134×10−60.34595×10−60.020.020200.12494×10−50.13883×10−50.020200.12494×10−50.13883×10−50.030.030440.28201×10−50.31337×10−50.030440.28201×10−50.31337×10−50.040.040790.50294×10−50.55888×10−50.040790.50294×10−50.55888×10−50.050.051230.78834×10−50.87605×10−50.051230.78834×10−50.87605×10−50.060.061780.000010.000010.061780.000010.000010.070.072420.000020.000020.072420.000020.000020.080.083160.000020.000020.083160.000020.000020.090.093990.000020.000030.093990.000030.000030.10.104930.000030.000040.104930.000030.00004Table 6
Differences between 9-iteration NIM with RK4 solutions in the sinusoidal input model.
∆=∣RK40.001−NIM∣t∆z1t∆z2t∆z3t21.3e-1402e-1649e-152e-165e-1661.0e-142e-168.4e-1581.85e-132.77e-141.690e-13102.33e-123.479e-131.974e-12121.726e-112.572e-121.462e-11149.251e-111.378e-117.863e-11163.933e-105.847e-113.349e-10181.404e-092.079e-101.196e-09204.366e-096.444e-103.722e-09Figure 6
Phase portraits using 9-iteration NIM for sinusoidal input.Figure 7
Comparison among 9-iteration NIM (solid line) and RK4 (circle) for (a)z1t, (b) z2t, and (c) z3tin the sinusoidal input model.
(a)(b)(c)All of the above discussions indicate that NIM is very accurate, requires less computation work, very easy to implement, and consumes less time to obtain good approximate solution in a wide interval for lake pollution problem.


## 6. Conclusions

In this study, the new iterative method (NIM) has been successfully implemented to gain the approximate analytic as well as numeric solutions of the lake pollution model and results compared with the traditional fourth-order Runge Kutta method (RK4) and ADM. The potential advantages of NIM are summarized below:(a)
This comparison shows that the solution procedure of the NIM is simpler than of ADM, LPDTM without any special type of polynomials or multipliers. Simple solution procedure, high accuracy and efficiency in computation are the advantages of the proposed method(b)
The accuracy of graphical and numerical results of NIM is very high in large time domain(c)
There are difficulties arising in calculating Adomian polynomials to handle the non-linear terms in ADM. It is necessary to calculate Lagrange multiplier in VIM in which the terms of the sequence become complex after several iterations; thus, finding more iterations in VIM becomes more laborious. Also, there is a need to construct a homotopy parameter in HPM and solving the corresponding algebraic equations will be somehow difficult. The advantages of NIM over other iterative methods (such as ADM, VIM, HPM) are its very straightforward in computer programming with less computation works and at the same time it provides fast convergent analytic-numeric solutions without any transformation, linearization, and discretizationThe outcomes introduced in this paper reveal that NIM is promptly pertinent to systems of ODEs including more complex physical models.



---

*Source: 6664307-2021-06-15.xml*