Integración numérica con MATLAB
Integrales simples
La función
>> f=@(x) -x.^2+14*x+21; >> q=integral(f,0,10) q = 576.6667
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
MATLAB dispone de la función
Creamos un script, en el que definimos la función anónima f(x,y) y llamamos a
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.
Vamos a 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
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.
Utilizamos la función MATLAB
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
res=integral(@g,0,2); fprintf('El valor de la integral es: %3.3f\n',res)
En la ventana de comandos corremos el script
El valor de la integral es: 961.181
Ejemplos
Comprobar los resultados de las integrales siguientes:
Comprobamos la primera integral. Cargamos el fichero
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
res=integral(@g1,0.1,0.5); 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: 0.0333
Comprobamos la segunda integral. Cargamos el fichero
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
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
Ejemplos en el curso de Física
Flexión de una viga en voladizo (II)
Pandeo de una barra delgada empotrada en un extremo
Fuerza de rozamiento proporcional al cuadrado de la velocidad
Vaciado de un depósito cerrado
Presión producida por la curvatura de una superficie
Función de distribución de Weibull