Metodo de Gauss-Legendre en MatLab

         El método Gauss - Lengendre, es una regla  o cuadratura perteneciente a las formulas de Newton- Cotes para la integración numérica con limites cerrados. 

          El código implementado a continuación se divide en dos funciones:

• La función “gaussleg”, esta función es la que recibe los parámetros con los cuales se van hacer los cálculos para generar la integración numérica.
• La función “gauss”, es la función que genera las abscisas y pesos.

Función “Gaussleg”:  

function [anss,abs1,wgt1] = gaussleg(f,a,b,gp,varargin)

format long

graf=(f);

f = fcnchk(f,’vectorized’);                

if nargin < 4, gp = 106; end              % Default numbero de puntos Gauss.

[abs1, wgt1] = Gauss(gp);  

Alphas = ((b-a)./2).*(abs1+1)+a;                     % Mapeo f a [-1,1].

anss = sum(wgt1.*f(Alphas))*(b-a)/2;

if isnan(anss) | isinf(anss) 

   fprintf([‘\n\t Incrementando el numbero de puntos Gauss ‘…

            por 1 porque esta a NaN.\n’])

   [anss,abs1,wgt1] = gaussleg(f,a,b,gp+1);

end

h= (b-a) /gp;

x= a+ (0:gp)*h;

y=eval(graf);

plot(x,y,’b-‘);

 grid;

  title(‘Gauss Legendre’),

  xlabel(‘x’);

  ylabel(‘f(x)’);

Funcion Gauss: 

function [x, w] = Gauss(n)

% Genera las abscisas y pesos para la  Cuadratura Gauss-Legendre.

x = zeros(n,1);                                            

w = x;

m = (n+1)/2;

for ii=1:m

    z = cos(pi*(ii-.25)/(n+.5));                        % estimado Inicial.

    z1 = z+1;

while abs(z-z1)>eps

    p1 = 1;

    p2 = 0;

    for jj = 1:n

        p3 = p2;

        p2 = p1;

        p1 = ((2*jj-1)*z*p2-(jj-1)*p3)/jj;       % El polinomial. Legendre 

    end

    pp = n*(z*p1-p2)/(z^2-1);                        % La L.P. Derivada.

    z1 = z;

    z = z1-p1/pp;

end

    x(ii) = -z;                                   % Construye las abscissas.

    x(n+1-ii) = z;

    w(ii) = 2/((1-z^2)*(pp^2));                     % Construye  los pesos.

    w(n+1-ii) = w(ii);

end