New handy and accurate approximation for the Gaussian integrals with applications to science and engineering

1 Introduction

Gaussian functions are considered among the most important special functions for science and engineering because they have a broad scope of application. Error function has applications in: chemical engineering, transfer phenomena [1], Newtonian fluids analysis to express analytical solutions of differential equations that model the flow near a wall suddenly set in motion and unsteady heat conduction in solids, microelectronics, digital and analogue signal processing, noise analysis for the transmission and communication signal protocols such as Phase Shift Keying (PSK) modulation [2, 3, 4, 5], among others.

Another Gaussian function is the normal distribution that has applications in the field of statistics and probability. This function allows to model several phenomena such as biological [6, 7, 8, 9], social [10, 11], psychological [12, 13], financial [14, 15, 16], in science and engineering [3, 4, 5, 17, 18, 19, 20, 21], among others.

Gaussian integrals can not be solved analytically employing traditional methods [22, 23]. However, there are some proposed methodologies in literature that allow to construct approximations. Therefore, in [24] are presented some approximative methods that are employed to evaluate Gaussian functions, such as Power Series Expansions, Rational Approximations, Continued Fraction Expansions, Approximation by Burr Distributions, Taylor series, among others. However, in order to obtain a good approximation such methods require to calculate higher order extra terms using an iterative procedure. This characteristics is a draw-back because the process is too cumbersome to be implemented by hand. Instead such approximative methods are programmed using specific languages such as Fortran, C++, among others. Another alternative to evaluate such integrals is by means of numerical integration, like Simpson’s or trapezoidal rules [25].

There exist other proposals to approximate Gaussian functions. For instance, in [24] a normal distribution integral analytical approximation is reported employing hyperbolic tangent function. Additionally, [26, 27] presents approximations related to error function expressed in terms of hyperbolic tangent and arctangent functions, while normal distribution integral was approximated with exponential and inverse tangent functions.

We propose to applicate the novel power series extender method (PSEM) [28, 29, 30] for the approximation of the three gaussian integrals, because PSEM exhibits a large domain of convergence and it does not requires a perturbation parameter [31, 32] or to calculate integrals [33, 34, 35] or the solution differential equations [36, 37, 38], as other aproximative methods.

In this work will assume that exact solution for Gauss functions are numerically solved employing Maple 15 and GNU Octave; in addition, one of the metrics that we will use to determine the precision of the proposed approximations will be estimated calculating the root-mean-squared (RMS) error, defined for our application as

where a and b are integration limits, E(t) is the relative error. The process of integration will be performed numerically utilizing Simpson’s rule 38 [39].

The organization of this article is as follows: Section 2 presents an introduction to Gauss integrals. Next, in Section 3 Power Series Extender Method is explained. Later, in Section 4 Gauss integrals approximations employing PSEM are presented. Three interesting applications to science and engineering are presented in Section 5. Discussion on the results for the approximations obtaining employing Gaussian integrals is presented in Section 6. Finally Section 7 shows the conclusions about this work.

2 A basic introduction to Gaussian integrals

This section presents some basis for the Gaussian distribution integral, the error function and the Cumulative distribution function (CDF).

2.1 The Gaussian distribution integral

The Gaussian integral, known as the Euler-Poisson integral, is Gaussian function exp(–π t²), which can not be analytically integrated [22, 23]

f(x)=∫0xexp⁡(−πt2)dt. (2)

The Gaussian distribution integral can be applied in: quantum mechanics to find the probability density for the fundamental state on the harmonic oscillator, the path integral formulation and the propagator for the harmonic oscillator. Figure 1 shows the well known Gaussian integral curve.

Figure 1

Gaussian integral (2).

2.2 Error function

Error function is classified as a special function and it is applied in the field of probabilistic, statistics, partial differential equations solutions, robotics, among others. Error function is defined as

erf(x)=2π∫0xexp⁡(−t2)dt, (3)

while complementary error function, erfc(x), is defined from error function as

erfc(x)=1−erf(x). (4)

Figure 2 presents the behaviour of (3).

Figure 2

Error function (3)

2.3 Cumulative distribution function

Normal distribution, also known as Gauss-Laplace distribution, is usually employed in statistics since aleatory processes tend to have similar behavior as this distribution. Several areas of knowledge present continuous aleatory variables that have a normal density distribution and their behavior exhibits a bell shape. It is noteworthy to mention the importance of normal distribution because many natural phenomena variables follow this model. In this tenor, in statistics, central limit theorem shows that under certain conditions (independent and identically distributed with finite variance), the sum of a large number of aleatory variables is approximately distributed as a normal function [3, 6, 7, 16, 21]. Normal probability density function is

p(x)=12πσexp⁡(−(x−μ)22σ2), (5)

where μ is the mean, σ the standard deviation and σ² the variance. Figure 3, shows that the graph is divided in function of standard deviation. For instance, –σ ≤ x ≤ σ interval has 66.36% of the area, which, in statistics is equal to have 66.36% of the data under study. In a similar fashion, –2σ ≤ x ≤ 2σ encloses 95.45% of the information and –3σ ≤ x≤ 3σ 99.73%. The probability is obtained by integrating from –∞ to x, which represents the area under the curve (see Figure 3). Thus, CDF is

Figure 3

Normal density function representation.

P(x)=12πσ∫−∞xexp⁡(−(x−μ)22σ2). (6)

Cumulative distribution function (CDF) is presented in Figure 4 for a particular case of σ and μ.

Figure 4

CDF normal (6), with σ=1 and μ=0.

3 Basic concept of PSEM method

In broad sense a nonlinear differential equation can be expressed as

having as boundary condition

where L and N are a linear operator and a nonlinear operator respectively; f(x), is a known analytic function; B, is a boundary operator; Γ, is the boundary of domain Ω; ∂u∂η denotes differentiation along the normal drawn outwards from Ω [40]. Next, we express the solution of (7) as a power series

where v_k (k = 0.1,2,…) are the coefficients of the power series.

It is important to notice that (9) can be obtained by some approximative method from literature HPM [36, 37, 38], VIM [34], TSM [41, 42], Taylor series, among others. Now, [29, 43] proposed that the solution for (7) can be written as a finite sum of functions in the general from

or

where u_i are constants to be determined by PSEM, f_i(x,u_i) are arbitrary functions, and n and 2n are the orders of approximations (10)1 and (11), respectively. We will denominate (10) and (11) as a trial function (TF).

Next, we calculate the Taylor series of (10) or (11), resulting in the power series:

respectively, where Taylor coefficients P_k are expressed in terms of parameters u_i.

Finally, we equate/match the coefficients of power series (12), (13) with (9) to obtain the values of u_i and substitute them into (10) or (11) to obtain the PSEM approximation. The proposed arbitrary functions can be functions, quotients of polynomials, transcendental functions, composite, products of transcendental functions, among others [29]. It is important to notice that PSEM convergence greatly depends on the proper selection of the trial function. Then, its necessary that the proposed TF can potentially describes the qualitative behavior of the solution of the nonlinear problem.

4 Approximation for the Gaussian integrals

This section presents approximation for the Gaussian integrals.

4.1 Approximation for the Gaussian distribution integral

Taylor series for the Gaussian Distribution Integral is

f(x)=x−13πx3+110π2x5−142π3x7+1216π4x9−⋯. (14)

We propose as TF to approximate the Gaussian Distribution Integral the following modified logistic function

f1~(x)=11+exp⁡(−c5x9−c4x7−c3x5−c2x3−c1x)−12. (15)

Next, we obtain the Taylor series of (15), resulting

f(x)=0+14c1x+(14c2−148c13)x3+(14c3+18c12c2+1480c15−18(c1c2+124c14)c1+14(−14c2+148c13)c12)x5+⋯. (16)

Equating coefficients from the respective x-powers of (14) and (16), the next nonlinear equation system results

x:14c1=1x3:14c2−148c13=−13πx5:14c3+18c12c2+1480c15−18(c1c2+124c14)c1+14(−14c2+148c13)c12=110π2x7:⋯=−142π3x9:⋯=1216π4 (17)

Solving the system (17) and substituting the values of c₁, c₂, c₃, c₄ y c₅ in (15), results

f1~(x)=11+exp⁡(ψ(x))−12,−∞<x<∞, (18)

where ψ(x) is the argument given by

ψ(x)=−(154π4−46642835π3+92845π2−2563π+10249)x9−(−221π3+15245π2−643π+2567)x7−(25π2−163π+645)x5−(−43π+163)x3−4x. (19)

In the literature, there are other expressions that allow to approximately calculate Gaussian function integrals. In [44] was proposed the following approximation

f2~(x)=12tanh⁡(39x2−1112arctan⁡(35x111)),−∞<x<∞. (20)

On one hand, Figure 5a shows a comparison among the exact solution and approximations (18) and (20). On the other hand, Figure 5b presents a comparison of relative error for both approximations resulting a notable lowest error for our proposal. In fact, RMS error of (18) for the interval [0,3] is 1.208×10^–5, while (20) shows a value of 1.742×10^–4; therefore, its error of our proposal is 14.42691943 times lower.

Figure 5

Gaussian distribution integral approximations and relative error.

4.2 Approximation for error function

Taylor series for error function is

erf(x)=2xπ−2x33π+x55π−x721π+x9108π−⋯. (21)

Next, the proposed trial function is the following modified logistic function

erf~1(x)=21+exp⁡(−c5x9−c4x7−c3x5−c2x3−c1x)−1. (22)

Repeating the PSEM procedure employed in the last section, we obtain the following approximation

erf~1(x)=21+exp⁡(ξ(x))−1,−∞<x<∞, (23)

where

ξ(x)=−15670(105π4−9328π3+116928π2−483840π+645120π9/2)x9 +2315(15π3−532π2+3360π−5760π7/2)x7−215(3π2−40π+96π5/2)x5+43(π−4π3/2)x3−4xπ1/2. (24)

In the literature there are several expressions to approximate the error function. For example, in [24, 27, 44] the following approximations were reported

erf~2(x)=tanh⁡(2x(1+0.089430x2)π),−∞<x<∞, (25)

erf~3(x)=1−exp⁡((−4π+0.14x2)x21+0.14x2),0≤x, (26)

erf~4(x)=tanh⁡(39x2π−1112arctan⁡(35x111π)),−∞<x<∞, (27)

respectively.

Figure 6a shows a comparison of exact solution and the above mentioned approximations. In addition, Figure 6b presents a comparison of relative error resulting a notable lower error for our proposal. Finally, Figure 6c depicts the relative error of our proposal in order to appreciate the error. Furthermore, calculating the RMS error for the interval [0,4] results: 1.392×10^–5, 2.248×10^–4, 1.551×10^–4, and 6.696×10^–5 for (23), (25), (26), and (27), respectively. It is important to notice that the error of our proposal (23) is several times lower than the other approximations.

Figure 6

Error function approximations and relative error.

4.3 Approximation for cumulative distribution function

The Taylor series expansion for CDF is

P(x)=12+2x2π−2x312π+2x580π−2x7672π+⋯. (28)

We propose the following TF function

P~1(x)=11+exp⁡(−c5x9−c4x7−c3x5−c2x3−c1x). (29)

After the application of PSEM procedure, we obtain

P~1(x)=11+exp⁡(ζ(x)),−∞<x<∞, (30)

where the argument ζ(x) is

ζ(x)=−2181440(105π4−9328π3+116928π2−483840π+645120π9/2)x9+22520(15π3−532π2+3360π−5760π7/2)x7−260(3π2−40π+96π5/2)x5 +23(π−4π3/2)x3−22π1/2x. (31)

Additionally, another way to approximate the CDF is obtained by employing the next TF

P~2(x)=c1+12tanh⁡(114arctan⁡(251250arctan⁡(c2x)+c3x)). (32)

Following PSEM presented in this article and solving for c₁, c₂ and c₃ yields to

P~2(x)=12+12tanh⁡(114arctan⁡(251250arctan⁡(−0.09839818659x)+0.1057907667x)),−∞<x<∞. (33)

Furthermore, in [44] was proposed the next approximation

P~3(x)=1exp⁡(−358x23+111arctan⁡(37x294))+1,−∞<x<∞. (34)

Figure 7a shows a comparison of the proposed PSEM approximation and the other ones reported, resulting a good agreement with the exact solution. However, Figure 7b can be observed the lowest relative error for our proposal. See Figure 7c for a detailed view of the relative error of (30). In addition, we calculated the RMS error over the interval [0,6] for the approximations, resulting: 6.756×10^–6, 3.509×10^–5, and 3.496×10^–5, for(30), (33) and (34), respectively. It clearly results that our proposal exhibits the lowest of all RMS errors.

Figure 7

Comparison of CDF and its approximations.

5 Applications to science and engineering

This section presents three case studies. The first case study belongs to the area of transport phenomena in chemical engineering [1]; the second one lies in the area of digital communications [2] in electronics engineering; and finally we propose an approximation for the incomplete gamma function.

5.1 Case 1: Error function in transport phenomena

A semi-infinite body in a liquid with constant density and viscosity limited under a solid plane xz horizontal. Initially the fluid is in a steady state; at t = 0, the solid surface begins to move in x positive direction with v_o velocity, as depicted in Figure 8a. The goal is to find the velocity v_x as a function of y and t. In the direction of x there is no pressure gradient neither gravity force and a laminar flow is considered.

Figure 8

Error function applied on time-dependent flow of Newtonian fluids.

In [1] the equation that models this system and the components of velocity at χ, y, z are v_χ = v_(y,t), v_y = 0, v_z = 0 is given by

∂vχ∂t=ν∂2vχ∂y2, (35)

where ν=μρ, with initial (IC) and boundary (BC) conditions

I.C.at t≥0,vχ=0,∀y,B.C.1at y=0,vχ=v0,∀t>0,B.C.2at y=∞,vχ=0,∀t>0.

Introducing an additional velocity ϕ=vχvo, equation (35) is re-written as

∂ϕ∂t=ν∂2ϕ∂y2, (36)

where

ϕ=ϕ(η), (37)

and

η=y4 νt. (38)

In [1] was presented the solution for (36) in terms of the function erf(x), resulting

vχ(y,t)v0=1−erf(y4νt). (39)

From (23) and (39), we obtain

vχ~(y,t)v0=2−21+exp⁡(ξ(y4νt)). (40)

Figure 8b depicts a comparison for the distribution of the viscous flow velocity of a fluid near a wall suddenly set in motion among (39) and (40), resulting a good agreement.

5.2 Case 2: Error function in digital communications

Bandwidth efficiency in Phase Shift Keying (PSK) modulation scheme increments utilizing M-PSK modulation [2]. In digital phase modulation, the M signals employed in PSK are represented as

Sm(t)=Re[g(t)exp⁡(jθm)exp⁡(j2πfct)], (41)

where S_m (t) is the waveform used to transmit information through the communication channel; g(t), is the time dependent signal pulse shape and θm=2π(m−1)M, where M are the possible phases of the carrier that transmits the information. If M = 2 we have a binary PSK. In this way, to have a more efficient use of the bandwidth, each signal element should represent more than a bit [45]. Figure 9 displays, for the case of BPSK for M = 2 the Signal space diagram of Binary Phase Shift Keying (BPSK) scheme, where S₁ and S₂ symbols are represented by binary numbers 1 and 0.

Figure 9

Signal space diagram of BPSK.

In digital transmission, Bit Error Rate (BER) is defined as the number of bits with errors divided by the total number of transmitted, received or processed bits for a determined period. The probability for BER in Binary Phase Shift Keying is

P(x)=12erfc(EbNo), (42)

where EbNo, is the Signal to Noise Ratio (SNR); E_b, the energy in one bit and N, is Additive White Gaussian Noise (AWGN). As the bit transmitted signal energy E_b increases, for a specific noise spectral density N_o; the messages corresponding to 1 and 0 symbols become more separated and the probability of error decreases [46]; as is depicted in Figure 10. By using (23), (42) can be rewritten as

Figure 10

Bit error probability curve for PSK modulation.

P~=1−11+exp⁡(ξ(EbNo)). (43)

Figure 10 depicts the bit error probability for PSK modulation obtained with (43). The simulation from this figure was obtained by replacing the error function by our PSEM approximation in the code Matlab/Octave published in [47] (see Appendix 1).

5.3 Case 3: Incomplete Gamma function

This case study presents the approximation for incomplete Gamma function [48] represented as

γ(a,x)=∫0xta−1exp⁡(−t)dt. (44)

Considering a=12,32,52,72 and 92; then (44) is simplified, resulting

γ(12,x)=πerf(π), (45)

γ(32,x)=−12(π erf(π)exp⁡(x)+2x)exp⁡(−x), (46)

γ(52,x)=−14(4x32−3π erf(π)exp⁡(x)+6x)exp⁡(−x), (47)

γ(72,x)=−18(8x52+20x32−15π erf(π)exp⁡(x)+30x)exp⁡(−x), (48)

γ(92,x)=−116(16x72+56x52+140x32−105π erf(π)exp⁡(x)+210x)exp⁡(−x), (49)

respectively.

Generally, for a=k2, where k < 0 and odd; the incomplete Gamma function is expressed in terms of error function. Therefore, the results from this case study can be extended to other values of a. Figure 11 presents the comparison of the exact solution of (45)-(49) against the results of using approximation erf~1(π), depicting a notable good agreement.

Figure 11

Exact Incomplete Gamma function versus (45)-(49).

6 Discussion

For the first case study, we obtained an accurate approximated solution for the distribution of the viscous flow velocity of a fluid near a wall suddenly set in motion with a relative error depicted in Figure 12. In this figure we can observe the points P_A and P_B that show the location where the relative error reaches peak magnitudes of –2.826×10^–6 and 3.372×10^–5, respectively. After point P_B the relative error tends to zero.

Figure 12

Relative error of (40).

For second case study, Figure 13 shows the analysis for the absolute error for the error function approximation applied to BER from BPSK. It is important to highlight the high accuracy of the results. The simulation for the second case study was performed utilizing the software Octave 4.0 replacing the built-in Matlab/Octave function erf(x) presented in [47] by (23). This case study demonstrated that the proposed approximate function erf~1(x) can be applied to digital signal processing with a high precision results as shown in Figure 13. Thus, our presented simulation shows a good performance which is in good agreement with the exact result presented in [47].

Figure 13

Bit error probability curve for PSK modulation in (43).

Figure 14 presents the absolute error of the approximate incomplete Gamma function for different values of a. In table 1 the RMS error over the interval [0,10] for every approximation of incomplete Gamma function is presented. There can be seen an increment for the RMS as a increases its value; this is because the relative error that erf~1(x) exhibits is scaled by a numeric coefficient that increases as a increases, resulting an unavoidable increment of the error. In order to mitigate this issue, the order of the PSEM approximation of erf(x) must be increased.

Figure 14

Absolute error.

It is important to point out that the present work shows that it is possible to obtain highly accurate asymptotic solutions employing PSEM. Likewise, from literature we know that truncated Taylor series have, regularly, a poor local convergence; however, the employment of truncated Taylor series as a part of PSEM, empowers the convergence of the approximated solutions; such as the ones presented in this work, which covers whole domain of the independent variable. In addition, it is relevant to note, that during the approximation process performed in this article, no integration procedure was performed; such as it happens with other approximate methods like HAM [37], HPM [49], VIM [34], classical perturbation method [31], Picard [50], Adomian [51, 52], among others, that require the application of analytical integration to obtain the approximate solutions depending of the number of iterations needed. Unfortunately, as we know from literature, many integrals does not have known exact solution at all producing a failure during the iterative steps (iterations) of such approximative methods. Finally, given that PSEM does not require a perturbation parameter [31] and it does not generate secular terms [49]; it can be concluded that this methodology presents a high potential for application in all areas of science and engineering where it is required the approximate solution for integrals, nonlinear differential equation, and special functions. Therefore, given the flexibility, applicability and novelty of PSEM method, we are proposing this method to be known also as Taylor-Leal Method (TLM) for the approximation of nonlinear problems where the Taylor expansion of the exact solution can be obtained.

7 Conclusions

This work proposed highly accurate PSEM approximations for the Gaussian distribution integral, error function and CDF. In fact, as the independent variable increases to infinity the absolute error tends to zero. It is worth to note that the approximations proposed in this work are handy and easy computable avoiding the application of data tables or numerical algorithms. Furthermore, we compared the proposed approximations with other reported in the literature finding that our proposals presented the lowest RMS error. Finally, three applications to science and engineering allow us to conclude that our proposed approximations exhibit a notable potential for the solution of different practical problems.