Cálculo integral. Regla de Barrow

En esta práctica pretendemos aprender a calcular integrales definidas e indefinidas (cuando sea posible) utilizando el programa Mathematica. Como motivación principal para la integral definida consideraremos el cálculo del área bajo la gráfica de una función positiva. Dedicaremos parte de la práctica a "ilustrar" gráficamente cómo la integral de Riemann es, desde su definición, una herramienta para calcular estas áreas y que la función de acumulación 
			F(x)=(∫_a^x) f(t)dt,
que aparece en el Teorema Fundamental del Cálculo Integral (T.F.C.I.) es realmente el "area acumulada" (con signo) de la función f (x). 
Aprenderemos los comandos que nos permiten calcular integrales definidas, calcular primitivas (con sus limitaciones) y concluiremos con la programación de los métodos de Trapecios y de Simpson para el cálculo aproximado de la integral definida.

En esta práctica aprenderemos los siguientes comandos:
         Integrate[ f [x], x]               	para el cálculo (simbólico) de la integral indefinida
         Integrate[ f [x], {x, a, b}]	para el cálculo (simbólico) de la integral definida (primitivas)
         NIntegrate[f [x], {x, a, b}] 	para el cálculo numérico de la integral definida

Integral definida. La integral de Riemann.

Introducción e interpretación gráfica

Recordamos rápidamente el contexto y la definición de la integral de Riemann (para una función continua y positiva en un intervalo cerrado y acotado). 
Dada   f : [a, b] -> R, una función continua y positiva pretendemos calcular el área que queda entre la gráfica de f, el eje OX y las rectas verticales x = a y x = b (para la ilustración hemos elegido f(x) = (x^3) - 4  (x^2) + 13 en el intervalo [0,4]).

Este es el origen y la definición de la integral definida (integral de Riemann). 
Para aproximar este área tomamos un natural n ϵ ℕ "muy grande" y partimos el intervalo en n trozos iguales usando los puntos de la forma (Subscript[x, i]) = i (( b - a)/n), con i=1,...,n. Al conjunto P={a = (Subscript[x, 0]) < (Subscript[x, 1]) < (Subscript[x, 2]) <  ...  < (Subscript[x, n]) = b} se le llama partición uniforme del intervalo [a,b]. Ahora se elige en cada subintervalo [(Subscript[x, i]), (Subscript[x, i + 1])] un punto cualquiera (Subscript[y, i]) y se aproxima el área de cada subintervalo por la del rectángulo de base (b-a)/n y  altura f((Subscript[y, i])). 
El área será más exacta cuantos más trocitos hagamos, y se alcanzará en el límite:
	Área = Integral de Riemann = Integral definida es:

			\(∫\_a\%b\( f(x)\) 
              d\ x\ \ \ \ \ \ \ \ es\ igual\ a\ \ \ \ \ \ \ \ \ \(\(\(\ \)\(lim\)\)\+\(n -> ∞\)\) ∑\+\(i = 1\)\%n\)f((Subscript[y, i]))*(((b - a))/n)

es decir, el límite de la suma de las áreas de los rectángulos conforme los intervalillos se van haciendo cada vez más pequeños. 

Ejecuta la celda siguiente para ilustrar gráficamente este límite.

<< 
      Graphics`FilledPlot`
\(integralMovie[
        func_, {x_, a_, b_}, {ymin_, ymax_}, sumRule_] := 
      Module[{apprx, f, stepfn, g, nrects, dx, plot}, 
\ \ \ f[xx_] := 
          func /. x -> xx; 
\ \ \ stepfn[xx_, n_] := 
          f[a + \((Floor[\(n \((xx - a)\)\)\/\(b - a\)] + 
                    sumRule)\) \(b - a\)\/n]; 
	\ g[xx_, n_] := 
          Max[{f[xx], stepfn[xx, n]}]; 
\ \ \ Do[nrects = 2\^k; 
          dx = \((b - a)\)/nrects // N; 
	
          apprx = Sum[
              f[a + \((i + sumRule)\) dx] dx, {i, 0, nrects - 1}]; 
	
          plot = FilledPlot[{f[x], stepfn[x, \ nrects], \ 
                g[x, \ nrects]}, {x, a, b + 0.001}, 
	
              MaxBend -> 0, Curves -> Front, 
	
              PlotStyle -> {{RGBColor[1, 0, 
                      0]}, {Thickness[ .005]}, {}}, 
              Fills\  -> {{{1, \ 2}, Hue[ .9,  .5, 1]}, {{2, Axis}, 
                    GrayLevel[ .85]}, \ {{1, 3}, 
                    RGBColor[ .7,  .8,  .9]}}, 
	
              PlotRange -> {ymin, ymax}, 
              DisplayFunction -> Identity]; 
          Show[plot, \ DisplayFunction -> $DisplayFunction], {k, 1, 8, 
            0.3}]]\)
integralMovie[
        (x^3) - 4  (x^2) - x + 13, {x, 0, 4}, {0, 13},  .5]

En esta animación hemos marcado en gris el área real, en rosa el error cometido por defecto y en azul el error por exceso; observamos como, al hacer la partición cada vez más fina, los errores van siendo cada vez menores.

Aproximando integrales definidas "a mano"

Basándonos simplemente en esta idea 
			\(∫\_a\%b\( f(x)\) 
              d\ x\ \ \ \ \ \ \ \ es\ igual\ a\ \ \ \ \ \ \ \ \ \(\(\(\ \)\(lim\)\)\+\(n -> ∞\)\) \(\(\ \)\(∑\)\)\+\(i = 1\)\%n\)f((Subscript[y, i]))*(((b - a))/n)
y usando Mathematica podríamos aproximar la integral definida de una función por la suma de la derecha con un valor de n suficientemente grande. 
 			\(TraditionalForm\`∫\_a\%b\( f(
            x)\)\ \ dx\ \ ~\ \(\(b - 
                  a\)\/n\) \((\ \(\(\(\ \)\(∑\)\)\+\(i = 1\)\%n\) \(f(
                x\_i)\))\)\)     (hemos elegido simplemente (Subscript[y, i])=(Subscript[x, i]))
A esta aproximación se le conoce como la Fórmula de los rectángulos. Apliquémosla a un ejemplo.
Elegimos la función f(x) = (Sqrt[1 - (x^2)]) para conocer de antemano el resultado de su integral en [-1,1], ya que corresponde a la mitad del área de un disco de radio 1, es decir (π/2). Si multiplicamos por 2 el resultado ha de ser π.
Varía tú el valor de n para ver cómo se acerca cada vez más al valor exacto

n = 
        200;
a = -1. ;
b = 1;
f[x_] := 
        2  Sqrt[1 - x^2];
aprox = (b - a)/n*
          Sum[f[-1 + i (b - a)/n], {i, 1, n}]

Plot[f[x], {x, -1, 1}, AspectRatio ->  .5];

Integral definida con Mathematica

¿Cómo calcular esta integral definida usando Mathematica? Es evidente que usar la definición en cada caso es poco práctico; existe un comando que realiza este cálculo directamente: Integrate.  Este comando se usa de la siguiente forma:
                                  Integrate[ f [x], {x, a, b}]	
y calcula (de manera simbólica) la integral definida de la función f en el intervalo [a,b]. Veamos unos ejemplos

Integrate[x^3 - 4  x^2 + 13, {x, 0, 4}]

Integrate[4/(1 + x^2), {x, 0, 1}]

Integrate[Log[x^2], {x, 0, 2}]

Integrate[1/x, {x, 0, 1}]

Como vemos, incluso detecta los errores si intentamos integrar algo que no es posible integrar. 
Pero el hacer "internamente" el cómputo de manera simbólica puede ser costoso y en ocasiones es más aconsejable trabajar de manéra numérica. Para ello usaremos el comando NIntegrate; este se usa de manera análoga:
                                  NIntegrate[ f [x], {x, a, b}]	
que calcula (usando internamente algún método numérico) la integral definida de la función f en el intervalo [a,b].

NIntegrate[x^3 - 4  x^2 + 13, {x, 0, 4}]
NIntegrate[4/(1 + x^2), {x, 0, 1}]
NIntegrate[Log[x^2], {x, 0, 2}]
NIntegrate[1/x, {x, 0, 1}]

Como vemos, también detecta que 1/x no es integrable en [0,1] debido a que su "método interno" no converge. En realidad, como veremos en teoría, el área que deja esta última función 1/x por debajo entre 0 y 1 es infinita, y por eso no es pocible hacer el cálculo.

Teorema fundamental del cálculo integral y regla de Barrow.
Cálculo de primitivas (integral indefinida)

Recordamos las definiciones y resultados fundamentales:
Dada una función continua f(x) en un intervalo [a,b]:
	- decimos que F es una primitiva de f cuando (F^( |))(x) = f(x) para cada x en [a,b].
	- la función de acumulación: F(x)=(∫_a^(x ))f(
              t) dt  es una primitiva de f. (T.F.C.I.)
	- Si G es una primitiva cualquiera de f, entonces (∫_a^b) f(x) dx = 
        G(b) - G(a), (regla de Barrow).

Cálculo de Primitivas

La Regla de Barrow proporciona una herramienta para calcular integrales definidas ¡siempre que se conozca una primitiva! La realidad es que, en la práctica, el cálculo teórico de primitivas suele ser más costoso que el cálculo numérico de la integral definida directamente, sin pasar por la primitiva. 
Mathematica es capaz de calcular simbólicamente la primitiva de una función (en la mayoría de los casos) e incluso reconoce algunas de las funciones que (sabemos) no tienen una primitiva explícita. El comando es también Integrate, pero sin indicar ningún intervalo de integración, es decir: Integrate[ f(x), x] calcula una primitiva de f con respecto a la variable x.
Obsérvense particularmente las salidas de las tres últimas integrales.

Integrate[x^3 - 4  x^2 + 13, x]

Integrate[4/(1 + x^2), x]

Integrate[Log[x^2], x]

Integrate[1/x, x]

Integrate[E^(x^2), x]

Integrate[(Sqrt[1 + (Cos[x])^2]), x]

Integrate[Sin[x]/x, x]

Sobre la función de acumulación

Veámos gráficamente cómo es la función de acumulación  F(x)=(∫_a^(x ))f(
              t) dt. 
Representamos (en negro) la función f(x) = (x/6) sen(x) en el intervalo [0,4π] y su función de acumulación F(x) (en azul). En la evolución vemos como al barrer áreas positivas (en verde) F crece, mientras que las "áreas negativas" (en rojo) al acumularse van restando y hacen decrecer a F. Por ello a F se le llama también función de acumulacion de área. Ejecuta la siguente celda para visualizarlo

areaMovie[func_, {x_, a_, b_}, {ymin_, ymax_}, frames_] := 
      Module[{curve, shade, 
          nada = Graphics[{GrayLevel[1], 
                Point[{2  b, 0}]}], Φgraph, f}, 
        f[t_] = func /. x -> t;  
        Do[curve = Plot[f[t], {t, a, b}, DisplayFunction -> Identity];  
          shade = If[a < x, 
              FilledPlot[{Max[f[t], 0], Min[f[t], 0]}, {t, a, x}, 
                Fills  -> {{{1, Axis}, Hue[ .42]}, {{2, Axis}, 
                      Hue[0]}}, Curves -> Front, 
                DisplayFunction -> Identity], 
              nada]; Φgraph = 
            If[a < x, 
              Plot[NIntegrate[f[s], {s, a, t}, AccuracyGoal -> 3], {t, 
                  a, x}, PlotDivision -> 4, 
                PlotStyle -> {Hue[ .7], Thickness[ .005]}, 
                DisplayFunction -> Identity], nada]; 
          Show[curve, shade, Φgraph, 
            Graphics[{Hue[ .7], PointSize[ .015], 
                Point[{x, 
                    NIntegrate[f[s], {s, a, x}, AccuracyGoal -> 3]}]}],  
            Graphics[{GrayLevel[ .5], Line[{{x, ymin}, {x, ymax}}]}],  
            PlotRange -> {{a, b}, {ymin, ymax}}, 
            DisplayFunction -> $DisplayFunction], {x, a, 
            b, (b - a)/(frames - 1)}]]
areaMovie[x/6 Sin[x], {x, 0, 4  π}, {-3, 3}, 20]

Ejercicios

Propuestos en clase por el profesor.