Numerical Integration: Gaussian Methods in 1D
The Gaussian integration, known also as the method of gaussian quadrature, is a numerical approximation of a definite integral of a function
in a general interval
. We will apply this approach to the computation of the integral of a polynomial function (which, for sufficiently big n, will give an exact result) and of general functions (approximate result). Gauss quadrature
Gauss quadrature formula is defined for the reference interval
as a weighted sume of the values of
at appropiately chosen points: 
.We need to choose both points and weights
. The points
are the zeros of Legendre polinomials
. This is a family of polynomial functions represented in the following figure according to the degree (only the first ones). See more information at this Wikipedia link. For i =1,..,5
The values and weights for these points in the Gauss quadrature formulas are shown in the following table

For example, if we choose n=3 the Gauss approximation formula is

.
Example:
Now we show an example of the different steps needed to compute the Gauss aprimation value of an integral. As you will see, due to that
is a polynomial function, the result is exact when we choose the proper number of Gauss points. Function to integrate
As an example, consider the polynomial function
, defined in the reference interval [-1,1] First we defined the function as an inline function using Matlab
f=@(x) x.^4-3*x.^2+1; %defined as an inline function
Order of integration (number of Gauss points)
When you fix a number of points n, the Gauss integration formula is exact for polynomials of degree 2n-1. In our example the polynomial degree is 4, then if we choose only two points n=2 the integral is still not exact. Therefore, the proper choose is n=3.
% because our function is a polynomial of order 4 we have to choose n=3
switch (n) %switch executes only one piece of the code according to the value of n
w=[1,1]; pG=[-1/sqrt(3), 1/sqrt(3)];
w=[5/9, 8/9, 5/9]; pG=[-sqrt(3/5), 0, sqrt(3/5)];
error('No data are defined for this value');
fprintf('n=%d , pG = %f %f %f ',n,pG)
n=3 , pG = -0.774597 0.000000 0.774597
The final Gauss formula
You can use a sequential formula
sumat = sumat + w(i)*f(pG(i));
or equivalently, a more compact form of the same procedure
Check the error
Here, because it is a simple example, we can compute the primitive of the function
and solve the integral using Barrow's rule. Then, we compare the result with the numerical value coming from the Gauss formula. primitiveF=@(x) x.^5/5-x.^3+x; %just to check the error
barrowRule=primitiveF(1)-primitiveF(-1);
errorInt = abs(barrowRule - intFcompact)
Generalization to integration over a general interval [a,b]
When the integration interval is not [-1,1] you have to make a change of variables in the intergral in order to apply the previous formulas to the new interval.
This means that the Gauss points, initially defined in the interval
, must be transformed into points on the new integration interval 
Example:
Compute the integral of the same function, but now in the interval
. So we have to use the change of variables for integrals. pGnew = a+(b-a)*(1+pG)/2; %we have already the value of pG
intFcompactNew = (w*f(pGnew)')*(b-a)/2 %the weight values are the same
Exercise 1:
According to the table included above, you have to use the function [w,pt] = gaussValues1D(n), included at the end of this script, to see how the approximation of the integral is improved (you must safe the gaussValues1D function in a new file named gaussValues1D.m).
f=@(x) x.^6-x.^4-3*x.^2+1; %defined as an inline function
primitiveF=@(x) x.^7/7-x.^5/5-x.^3+x; %just to check the error
barrowRuleNew=primitiveF(3)-primitiveF(-2);
for n=1:5 %for n=6 must return an error
[w,pG]=gaussValues1D(n); %function to be implemented
pGnew = a+(b-a)*(1+pG)/2; %we have already the value of pG
intFcompactNew = (w*f(pGnew)')*(b-a)/2; %the weight values are the same
errorInt = abs(barrowRuleNew - intFcompactNew);
fprintf('for n= %d the error is =%e \n',n,errorInt)
end
for n= 1 the error is =2.446987e+02
for n= 2 the error is =1.769180e+02
for n= 3 the error is =2.790179e+01
for n= 4 the error is =2.842171e-13
for n= 5 the error is =4.263256e-13
Exercise 2:
Using the function gaussValues1D, approximate the value of the integral in
of the function
. Use different values of n and show the errors comparing with the true integral value 
Sol:
n, error
Function gaussValues 1D
function [w,pt]=gaussValues1D(n)
w=[1,1]; pt=[-1/sqrt(3), 1/sqrt(3)];
w=[5/9, 8/9, 5/9]; pt=[-sqrt(3/5), 0, sqrt(3/5)];
w=[(18+sqrt(30))/36, (18+sqrt(30))/36, (18-sqrt(30))/36, (18-sqrt(30))/36];
pt=[-sqrt(3/7-2/7*(sqrt(6/5))), sqrt(3/7-2/7*(sqrt(6/5))),...
-sqrt(3/7+2/7*(sqrt(6/5))), sqrt(3/7+2/7*(sqrt(6/5)))];
w=[(322+13*sqrt(70))/900, (322+13*sqrt(70))/900, 128/225,...
(322-13*sqrt(70))/900,(322-13*sqrt(70))/900];
pt=[-1/3*sqrt(5-2*sqrt(10/7)), 1/3*sqrt(5-2*sqrt(10/7)),0,...
-1/3*sqrt(5+2*sqrt(10/7)), 1/3*sqrt(5+2*sqrt(10/7))];
error('No data are defined for this value');
(c) Numerical Factory 2023