Unidad V: Integración numérica.

Unidad V: Integración numérica.

En análisis numérico, la integración numérica constituye una amplia gama de algoritmos para calcular el valor numérico de una integral definida y, por extensión, el término se usa a veces para describir algoritmos numéricos para resolver ecuaciones diferenciales. El término cuadratura numérica (a menudo abreviado a cuadratura) es más o menos sinónimo de integración numérica, especialmente si se aplica a integrales de una dimensión a pesar de que para el caso de dos o más dimensiones (integral múltiple) también se utiliza.

Hay varias razones para llevar a cabo la integración numérica. La principal puede ser la imposibilidad de realizar la integración de forma analítica. Es decir, integrales que requerirían de un gran conocimiento y manejo de matemática avanzada pueden ser resueltas de una manera más sencilla mediante métodos numéricos. Incluso existen funciones integrables pero cuya primitiva no puede ser calculada, siendo la integración numérica de vital importancia. La solución analítica de una integral nos arrojaría una solución exacta, mientras que la solución numérica nos daría una solución aproximada. El error de la aproximación, que depende del método que se utilice y de qué tan fino sea, puede llegar a ser tan pequeño que es posible obtener un resultado idéntico a la solución analítica en las primeras cifras decimales

En análisis numérico la regla del trapecio es un método de integración, es decir, un método para calcular aproximadamente el valor de una integral definida. La regla se basa en aproximar el valor de la integral de  por el de la función lineal, que pasa a través de los puntos  y . La integral de ésta es igual al área del trapecio bajo la gráfica de la función lineal.

Integrales dobles

Queremos calcular una integral doble de la función f(x,y) en la región rectangular de abscisas a y b y de ordenadas c y d, tal como se muestra en la figura

∫ba∫dcf(x,y)⋅dy⋅dx∫ab∫cdf(x,y)·dy·dx

MATLAB dispone de la función integral2 que realiza esta tarea. Supongamos que queremos calcular la integral doble

∫21∫30x2y⋅dy⋅dx∫12∫03x2y·dy·dx

Creamos un script, en el que definimos la función anónima f(x,y) y llamamos a integral2pasándole la función, los límites de integración de la variable x, y los límites de integración de la variable y.

f=@(x,y) x.^2.*y;

res=integral2(f,1,2,0,3);

fprintf('El valor de la integral es: %2.3f\n',res)

En la ventana de comandos corremos el script

 El valor de la integral es: 10.500

Como puede comprobarse fácilmente haciendo la integral doble.

Para calcular la integral doble de la función f(x,y) en la región no rectangular comprendida entre las curvas y=c(x), y=d(x) y las rectas x=a y x=b

b∫ad(x)∫c(x)f(x,y)dy⋅dx∫ab∫c(x)d(x)f(x,y)dy·dx

Al integrar la función f(x,y) respecto de la variable y entre los límites c(x) y d(x) nos queda la función g(x) que integramos entre los límites a y b.

∫ba(∫d(x)c(x)f(x,y)⋅dy)dx=∫bag(x)⋅dx∫ab(∫c(x)d(x)f(x,y)·dy)dx=∫abg(x)·dx

Utilizamos la función MATLAB quad, para resolver ambas integrales, tal como vemos en el siguiente ejemplo.

∫20∫x2+xx3−x(x3y4+xy2)⋅dy⋅dx=2∫0(∫x2+xx3−x(x3y4+xy2)⋅dy)dx∫02∫x3−xx2+x(x3y4+xy2)·dy·dx=∫02(∫x3−xx2+x(x3y4+xy2)·dy)dx

Definimos la función g(x) del siguiente modo

function z=g(x)

    n=length(x);

    z=zeros(size(x));

    for j=1:n

    %límites de y

        c=x(j)^3-x(j);

        d=x(j)^2+x(j);

    %integrando f(x,y)

        f=@(y) x(j)^3*y.^4+x(j)*y.^2;

        z(j)=integral(f,c,d);

    end

end

En el script llamamos a la función quad para resolver la integral de g(x)

res=integral(@g,0,2);

fprintf('El valor de la integral es: %3.3f\n',res)

En la ventana de comandos corremos el script integral_doble1

   El valor de la integral es: 961.181

Ejemplos

Comprobar los resultados de las integrales siguientes:

0.5∫0.1x2∫x3ey/xdy⋅dx≈0.03331∫02x∫x(x2+y3)dy⋅dx≈1.0001∫0.10.5∫x3x2ey/xdy·dx≈0.0333∫01∫x2x(x2+y3)dy·dx≈1.0001

Comprobamos la primera integral. Cargamos el fichero g.m en el editor y lo guardamos (File/Save As..) con nombre g1.m. Cambiamos el nombre de la función g1 en vez de g. Modificamos los límites c y d de la integral así como la definición de la función f(x,y)

function z=g1(x)

    n=length(x);

    z=zeros(size(x));

    for j=1:n

    %límites de y

        c=x(j)^3;

        d=x(j)^2;

    %integrando f(x,y)

        f=@(y) exp(y./x(j));

        z(j)=integral(f,c,d);

    end

end

Modificamos el script para llamar a la función integral para calcular la integral de g1(x) entre los nuevos límites de la integral a y b

res=quad(@g1,0.1,0.5);

fprintf('El valor de la integral es: %1.4f\n',res)

En la ventana de comandos corremos el script integral_doble2

El valor de la integral es: 0.0333

Comprobamos la segunda integral. Cargamos el fichero g1.m en el editor y lo guardamos (File/Save As..) con nombre g2.m. Cambiamos el nombre de la función g2 en vez de g1. Modificamos los límites c y d de la integral así como la definición de la función f(x,y)

function z=g1(x)

    n=length(x);

    z=zeros(size(x));

    for j=1:n

    %límites de y

        c=x(j);

        d=2*x(j);

    %integrando f(x,y)

        f=@(y) x(j)^2+y.^3;

        z(j)=integral(f,c,d);

    end

end

Modificamos el script para llamar a la función integral para calcular la integral de g2(x) entre los nuevos límites de la integral a y b

res=integral(@g2,0,1);

fprintf('El valor de la integral es: %1.4f\n',res)

En la ventana de comandos corremos el script

El valor de la integral es: 1.0000

comprendida entre las curvas y=c(x), y=d(x) y las rectas x=a y x=b