# New Modification on Heun’s Method Based on Contraharmonic Mean for Solving Initial Value Problems with High Efficiency

Authors: Abushet Hayalu Workie
Journal: Journal of Mathematics
 (2020)

---

## Abstract

In this paper, small modification on Improved Euler’s method (Heun’s method) is proposed to improve the efficiency so as to solve ordinary differential equations with initial condition by assuming the tangent slope as an average of the arithmetic mean and contra-harmonic mean. In order to validate the conclusion, the stability, consistency, and accuracy of the system were evaluated and numerical results were presented, and it was recognized that the proposed method is more stable, consistent, and accurate with high performance.


---

## Body

## 1. Introduction

Differential equations, either ordinary derivatives or partial derivatives, are equations which contain derivatives. An ordinary differential equation together with initial condition is an initial value problem (IVP) which specifies the value of the unknown function at a given point in the domain. There are a lot of physical problems in science and engineering which exist in the form of differential equations and are also commonly used in physics, chemistry, biology, and economics [1].To determine the solution of differential equations, there are different analytical methods available. In certain cases, however, analytical methods are not capable of solving such complicated or complex differential equations. Numerical methods are used to achieve the solution to the complicated differential equations [2, 3]. With the help of computer programming, numerical methods are very valuable tools for solving complex problems very quickly.Numerous numerical methods for solving ordinary differential equations (ODEs) with initial value problems (IVPs) have been developed by many researchers. Many authors have attempted to solve initial value problems (IVP) to obtain high accuracy rapidly by using a numerous methods, such as Euler’s method and Heun’s method, and some other methods. Euler’s method uses the line tangent to the function at the beginning of the interval as an estimate of the slope of the function over the interval. However, Heun’s method considers the tangent lines to the solution curve at both ends of the interval. And some of them have attempted to enhance these precision methods, where others have improved these methods for better accuracy, stability, and consistency [4–7]. Some improvements have been made from time to time in numerical methods to get better performance according to our needs.In this paper, Heun’s method and its modification are applied for solving ordinary differential equation in initial value problems. The numerical results are very encouraging. This demonstrates improved performance and better accuracy compared to other well-known methods of second order present in the literature.


## 2. Methodology of Research of the Proposed Method

Let us consider the first-order differential equation with initial value problem (IVP):(1)dydx=fx,ywith initial conditionyx0=y0.In order to solve (1), the most known and simplest method is the explicit Euler’s (Forward Euler) method with step length h and is given by(2)yn+1p=yn+hfxn,yn.This formula comes from approximating the derivativey0 at x=xn by a forward difference.Improved Euler’s method (Heun’s method) has a slope which is the average of the slopes of the tangents to the integral curve at the endpoints of the intervalxi,xi+1 given by(3)yn+1c=yn+h2fxn,yn+fxn+1,yn+1p.Improved Euler’s method is also known as Heun’s method. In this method, Euler’s method is used as predictor step and a mean slope is used as a corrector step.In [8–10], a second stage of the second-order harmonic mean method is given by(4)yn+1=yn+hK12+K22K1+K2whereK1=fxn,ynandK2=fxn+1,yn+1.Now, the new proposed method assumes the slope as the average of the arithmetic mean and second-order harmonic mean of the slopes of the tangents which is given by equation (5) as(5)yn+1c=yn+h2K1+K22+K12+K22K1+K2,where K1=fxn,ynandK2=fxn+1,yn+1p=fxn+h,yn+hfxn,yn.


## 3. Stability Analysis of the Proposed Method

The stability analysis of the proposed method can be obtained using Dahlquist’s test problem [11]:(6)dydx=λyn;yx0=y0,λ∈ℂ.Heun’s method and the new proposed (improved Heun’s) method require two evaluations offx,y per step, while Euler’s method requires only one. The local truncation error with the improved Euler method is Oh3, rather than Oh2 as with Euler’s method.IfK1=fxn,ynandK2=fxn+1,yn+1p, then the new proposed (improved Heun’s) method is the average of the artihmetic mean and the second order contra harmonic mean which is given as(7)yn+1c=yn+h2K1+K22+K12+K22K1+K2.In order to check the stability of the proposed method, substitute (6) into (7), and we have(8)yn+1c=yn+h2λyn+λyn+1p2+λyn2+λyn+1p2λyn+λyn+1p=yn+h2λyn+λyn+hλyn2+λyn2+λyn+hλyn2λyn+λyn+hλyn=yn+h22λyn+hλ2yn2+λyn+λ1+λh2yn2+hλ=1+λh2+λ2h24+λh+λh1+λh22+hλyn.If we letz=hλ, we have(9)Gz=yn+1yn=1+z2+z24+z+z1+z22+z=3+z2+54z2−4z+2;z≠−2.The stability region is(10)Gz=3+z2+54z2−4z+2≤1.Here, in Figure1, the shaded region is the region of unstability of the proposed method. In other words, the proposed method is stable in the unshaded region.Figure 1
Stability region of the new proposed method.


## 4. Consistency Analysis

In order to check the consistency of the initial value problem (IVP), the numerical formula expressed by [11], equation (1) can be written in the form of(11)yn+1=yn+hϕxn,yn,h,and will be consistent if limh⟶0ϕxn,yn,h=fxn,yn.In our case, in order to check the consistency of the proposed method of equation (7), we have(12)ϕxn,yn,h=12K1+K22+K12+K22K1+K2=12fxn,yn+fxn+1,yn+1p2+fxn,yn2+fxn+1,yn+1p2fxn,yn+fxn+1,yn+1p.From this,(13)limh⟶0ϕxn,yn,h=limh⟶012fxn,yn+fxn+1,yn+1p2+fxn,yn2+fxn+1,yn+1p2fxn,yn+fxn+1,yn+1p=fxn,yn.Therefore, the proposed method is consistent.


## 5. Local Truncation Error

The local truncation error (LTE) of a numerical method is an estimate of the error introduced in a single iteration of the method, assuming that everything fed into the method was perfectly accurate.The Taylor series expansion ofyxn+1 is(14)yxn+1=yxn+h=yxn+hy′xn+h22y〞xn+h33!y‴xn+Oh4,(15)yxn+1=yxn+hf+h22fx+ffy+h33!fxx+2ffxy+f2fyy+ffy2+fxfy+Oh4.The Taylor series expansions of the slopesK1 and K2 of the tangents are(16)K1=fxn,yn=f,K2=fxn+1,yn+1p=fxn+1,yn+hftn,yn=f+hfx+ffy+h22fxx+2ffxy+f2fyy+Oh3,respectively.And then, the sum of arithmetic mean and contra harmonic mean ofK1 and K2 is given as(17)K1+K22+K12+K22K1+K2=2f+hfx+ffy+h242f2fyy+4ffxy+2fxx+ffx2+2f2fxfy+f3fy2−1f2fx+ffy2+Oh3.Substituting (17) into the proposed equation (7), we obtain(18)yn+1=yn+h22f+hfx+ffy+h242f2fyy+4ffxy+2fxx+ffx2+2f2fxfy+f3fy2−1f2fx+ffy2+Oh3.After comparing equations (15) and (18), the local truncation error (LTE) is given by(19)LTE=yxn+1−yn+1and hence LTE=h3−112fxx−112f2fyy+13ffxy+16fxfy−18ffx2−14f2fxfy−18f3fy2+18f2fx+ffy2+Oh4.Therefore, the new proposed method has a local truncation error of orderOh3 and hence, the order of the proposed method is of order 2.


## 6. Numerical Results and Discussion

In this section, numerical solutionsyn of yx for the initial value problem (IVP) using Heun’s method and the new Proposed method has been discussed in Examples 1–3, and Tables 1–6 show the numerical solutions and errors made, while Figures 2–4 show the graphs of the numerical solutions using Heun’s method for different step lengths and the exact solution. Figures 5 and 6 show the graphs of the numerical solutions using the new proposed method, and Figures 7 and 8 show the errors made using the new proposed method and Heun’s method based on different h− values, and the results are computed with the help of MATLAB software [12].Table 1
Numerical solutions and exact solution for Example1.
Heun’s methodThe proposed methodExact solutionXh = 0.1h = 0.05h = 0.1h = 0.05Exact0.01.01.01.01.01.00.10.8200409370.8190505720.8189250180.8187883070.8187512210.20.6727344500.6710864550.6705490800.6706359420.6705881740.30.5525976430.5505438780.5498040500.5499774520.5499229800.40.4551606370.4528906160.4519982580.4522555450.4522046690.50.3766812510.3743357470.3733457810.3736688350.3736275570.60.3139709200.3116522390.3106018290.3109804370.3109529040.70.2642876110.2620676240.2609881680.2614336010.2613989470.80.2252677020.2231942810.2221124710.2225835250.2225707210.90.1948795010.1929817570.1919210280.1924279270.1924120381.00.1713880700.1696806730.1686614910.1691652900.169169104Table 2
Errors made by the stated numerical methods for Example1.
No. of iterationError made during Heun’s method whenError made during the proposed method whenh=0.1h=0.05h=0.1h=0.050000010.0012897162.99351e − 0041.73797e − 0043.7086e − 00520.0021462764.98281e − 0043.90940e − 0054.7768e − 00530.0026746636.20898e − 0041.18930e − 0045.4472e − 00540.0029559686.85947e − 0042.06411e − 0045.0876e − 00550.0030536947.08190e − 0042.81776e − 0044.1278e − 00560.0030180166.99335e − 0043.51075e − 0042.7533e − 00570.0028886646.68677e − 0044.10779e − 0043.4654e − 00580.0026969816.23560e − 0044.58250e − 0041.2804e − 00590.0024674635.69719e − 0044.91010e − 0041.5889e − 005100.0022189665.11569e − 0045.07613e − 0043.8140e − 006Table 3
Numerical solutions and exact solution for Example2.
Heun’s methodThe proposed methodExact solutionxh = 0.1h = 0.05h = 0.1h = 0.05Exact0.0−1.0000000−1.000000−1.0000000−1.0000000−1.00000000.1−0.9150000−0.914125−0.9140760−0.9144658−0.91451220.2−0.8570000−0.856392−0.8560341−0.8560892−0.85619220.3−0.8235250−0.822380−0.8222120−0.8224399−0.82245470.4−0.8122350−0.810974−0.8108372−0.8109677−0.81096010.5−0.82102288−0.819677−0.8191856−0.8195801−0.8195920Table 4
Errors made by the stated numerical methods for Example2.
Stage of iterationError made during Heun’s method whenError made during the proposed method whenh=0.1h=0.05h=0.1h=0.050000014.8780e−0043.8720e−0044.3620e−0044.64000e−00528.0780e−0041.9980e−0041.5820e−0041.23999e−00431.0703e−0031.3900e−0052.4269e−0041.47999e−00541.2749e−0031.3900e−0051.2280e−0047.60000e−00651.4308e−0038.5000e−0054.0640e−0041.18999e−005Table 5
Numerical solutions and exact solution for Example3.
Heun’s methodThe proposed methodExact solutionxh = 0.5h = 0.2h = 0.5h = 0.25Exact0.0111110.50.6050.58382495110.55180379750.56956461930.57216187271.00.44056250.43024902030.38988033000.41730092940.42035038451.50.5286750.52028910380.50600270000.51534904990.50915642062.01.1785436251.27442117591.20756012861.31591476891.3056051721Table 6
Errors made by the stated numerical methods for Example3.
Stage of iterationError made during Heun’s method whenError made during the proposed method whenh=0.5h=0.25h=0.5h=0.250000010.0328381273000.011663078400.02035807520.0025972534020.020212115500.0098986358000.030470054500.0030494551030.0195185794000.011132683200.003153720600.00619262930040.12706154710.03118399620.098045043500.01030959680Figure 2
Numerical solution using Heun’s method and exact solution for Example1.Figure 3
Numerical solutions using Heun’s method and exact solution for Example2.Figure 4
Numerical solutions and exact solution for Example3.Figure 5
Numerical solution using the proposed method and exact solution for Example1.Figure 6
Numerical solutions using the proposed method and exact solution for Example2.Figure 7
Errors made by the stated methods for Example1.Figure 8
Errors made by the stated methods for Example3.Example 1.
Solve the initial value problem (IVP)y′=x3e−2x−2y;y0=1 using Heun’s (improved Euler) method and the proposed method with step sizes h=0.1 and h=0.05 and compare the errors.Example 2.
Solvey′=−2x−y;y0=−1 using Heun’s (improved Euler) method and the proposed method for h=0.1andh=0.05.Example 3.
Solve the IVPy′=yx2−1.2y;y0=1 for h=0.5 and h=0.25.
Solution: for Example1, Table 1 shows results using Heun’s method and the new proposed method with step sizes h=0.1 and h=0.05 and the exact solution, while Table 2 shows the error made using Heun’s method and the new proposed method to find approximate values of the solution of the initial value problem.
For Example2, Table 3 shows results using Heun’s method and the new proposed method with step sizes h=0.1 and h=0.05 and the exact solution, while Table 4 shows the error made using Heun’s method and the new proposed method to find approximate values of the solution of the initial value problem.
For Example3, Table 5 shows results using Heun’s method and the new proposed method with step sizes h=0.5 and h=0.25 and the exact solution, while Table 6 shows the error made using Heun’s method and the new proposed method to find approximate values of the solution of the initial value problem.


## 7. Conclusion

This research paper shows the possibility of constructing a new proposed method on Heun’s method based on second-order contra-harmonic mean to solve ordinary differential equations (ODEs) with initial value problems (IVPs). The result of the new proposed method has been compared with Modified Euler method and Heun’s method (Improved Euler method). Tables1 and 3 demonstrate the numerical solutions of Heun’s method and the new proposed method with step lengths h=0.1 and h=0.05, and the results show that the new proposed method at h=0.1 is by far better than Heun’s method at h=0.05. Tables 2 and 4 show the absolute errors of the numerical solutions obtained under Tables 1 and 3, respectively. Table 5 demonstrates the numerical solution of Heun’s method and the new proposed method with step lengths h=0.5 and h=0.25, and Table 6 shows the absolute errors of the numerical solutions obtained under Table 5.The new proposed method has been constructed or formed by combining the slopesK1 and K2 and using the slope as the average of the arithmetic mean and contra-harmonic mean. And also the stability and consistency have been analyzed, and the local truncation error was computed and found that the new proposed method is stable, consistent, and second-order accurate.



---

*Source: 6650855-2020-12-07.xml*